Modèle linéaire généralisé et Choix de modèles
Philippe Real
30/8/2019
#### Instalation des packages.
#for manipulate data (transform to dataframe)
install.packages("bindrcpp")
install.packages("pillar")
install.packages("tidyverse")
install.packages("tibble")
install.packages("sm")
install.packages("leaps")
install.packages("stringr", dependencies=TRUE)
install.packages("readxl")
install.packages("VIM")
install.packages("FactoMineR")
install.packages("ggplot2")
install.packages("factoextra")
install.packages("corrplot")
install.packages("devtools")
#install.packages("usethis")
install.packages("plyr")
install.packages("rgl")
install.packages("rJava")
install.packages("glmulti")
install.packages("questionr")
install.packages("faraway")Introduction
Après avoir chargé les données au §1 on construit le modèle saturé au §2, modèle : m_sature. Ensuite on effectue une analyse PCA au §3 pour mieux appréhender les données. De cette analyse des données on va construire un modèle alternatif au modèle saturé, le modèle PCA: m_PCA Ensuite au §4 vient la sélection des variables. Pour celà on va utiliser le package MASS et la méthode step. Tout d’abord à partir du modèle saturé: m_sature et aussi à partir du modèle m_PCA, comme alternative. La pré-selection effectuée avec l’analyse PCA va aussi nous permettre d’utiliser d’autres méthodes. Comme la méthode leaps du package leaps. En effet bon nombre de méthodes ne sont applicables qu’à partir de 30 variables. Puis au §5 on essaiera une analyse croisée et comparaison de modèles. Après avoir choisi et validé notre modèle au §6, on effectuera la prédiction au §7. Ceci conduira à la production, à partir du modèle retenu, d’un fichier de prédiction avec nos données de test.
1. Lecture et préparation des données, premières analyses
Le nom des colonnes est trop long gêne dans les graphiques et la lecture des données. On va renommer les colonnes :
| Nouveau Nom | Ancien Nom |
|---|---|
| CC.high_max | High.Cloud.Cover.daily.max..high.cld.lay. |
| CC.high_mean | High.Cloud.Cover.daily.mean..high.cld.lay. |
| CC.high_min | High.Cloud.Cover.daily.min..high.cld.lay. |
| CC.low_max | Low.Cloud.Cover.daily.max..low.cld.lay. |
| CC.low_mean | Low.Cloud.Cover.daily.mean..low.cld.lay. |
| CC.low_min | Low.Cloud.Cover.daily.min..low.cld.lay. |
| MSL_max | Mean.Sea.Level.Pressure.daily.max..MSL. |
| MSL_mean | Mean.Sea.Level.Pressure.daily.mean..MSL. |
| MSL_min | Mean.Sea.Level.Pressure.daily.min..MSL. |
| CC.mid_max | Medium.Cloud.Cover.daily.max..mid.cld.lay. |
| CC.mid_mean | Medium.Cloud.Cover.daily.mean..mid.cld.lay. |
| CC.mid_min | Medium.Cloud.Cover.daily.min..mid.cld.lay. |
| Humidity_max | Relative.Humidity.daily.max..2.m.above.gnd. |
| Humidity_mean | Relative.Humidity.daily.mean..2.m.above.gnd. |
| Humidity_min | Relative.Humidity.daily.min..2.m.above.gnd. |
| Rad_sum | Shortwave.Radiation.daily.sum..sfc. |
| Snow_sum | Snowfall.amount.raw.daily.sum..sfc. |
| Sunshine_sum | Sunshine.Duration.daily.sum..sfc. |
| Temp_max | Temperature.daily.max..2.m.above.gnd. |
| Temp_mean | Temperature.daily.mean..2.m.above.gnd. |
| Temp_min | Temperature.daily.min..2.m.above.gnd. |
| Tot.CC_max | Total.Cloud.Cover.daily.max..sfc. |
| Tot.CC_mean | Total.Cloud.Cover.daily.mean..sfc. |
| Tot.CC_min | Total.Cloud.Cover.daily.min..sfc. |
| Tot.Prcp_sum | Total.Precipitation.daily.sum..sfc. |
| WDir.10m_mean | Wind.Direction.daily.mean..10.m.above.gnd. |
| WDir.80m_mean | Wind.Direction.daily.mean..80.m.above.gnd. |
| WDir.900mb_mean | Wind.Direction.daily.mean..900.mb. |
| WGust_max | Wind.Gust.daily.max..sfc. |
| WGust_mean | Wind.Gust.daily.mean..sfc. |
| WGust_min | Wind.Gust.daily.min..sfc. |
| WS.10m_max | Wind.Speed.daily.max..10.m.above.gnd. |
| WS.80m_max | Wind.Speed.daily.max..80.m.above.gnd. |
| WS.900mb_max | Wind.Speed.daily.max..900.mb. |
| WS.10m_mean | Wind.Speed.daily.mean..10.m.above.gnd. |
| WS.80m_mean | Wind.Speed.daily.mean..80.m.above.gnd. |
| WS.900mb_mean | Wind.Speed.daily.mean..900.mb. |
| WS.10m_min | Wind.Speed.daily.min..10.m.above.gnd. |
| WS.80m_min | Wind.Speed.daily.min..80.m.above.gnd. |
| WS.900mb_min | Wind.Speed.daily.min..900.mb. |
Il y a certaienement une manière plus “intelligente”, par exemple à l’aide de %>% des filter et des recherche de caractères… En tout cas faire très attention au rennomage pour garder la cohérence des variables et des données. A bien vérifier donc. Ce choix est certainement discutable mais les noms de variables sont vraiment trop longs.
meteoTrainAll<-renameCol(meteoTrainAll)
summary(meteoTrainAll)## X Year Month Day
## Min. : 1.0 Min. :2010 Min. : 1.000 Min. : 1.00
## 1st Qu.: 700.5 1st Qu.:2012 1st Qu.: 4.000 1st Qu.: 8.00
## Median :1438.0 Median :2014 Median : 6.000 Median :16.00
## Mean :1446.9 Mean :2014 Mean : 6.493 Mean :15.73
## 3rd Qu.:2183.5 3rd Qu.:2016 3rd Qu.: 9.000 3rd Qu.:23.00
## Max. :2941.0 Max. :2018 Max. :12.000 Max. :31.00
## Hour Minute Temp_mean Humidity_mean MSL_mean
## Min. :0 Min. :0 Min. :-7.100 Min. :41.12 Min. : 978.5
## 1st Qu.:0 1st Qu.:0 1st Qu.: 6.695 1st Qu.:65.11 1st Qu.:1012.6
## Median :0 Median :0 Median :12.575 Median :71.88 Median :1017.0
## Mean :0 Mean :0 Mean :12.227 Mean :71.45 Mean :1017.0
## 3rd Qu.:0 3rd Qu.:0 3rd Qu.:17.640 3rd Qu.:78.42 3rd Qu.:1021.6
## Max. :0 Max. :0 Max. :29.590 Max. :95.25 Max. :1038.6
## Tot.Prcp_sum Snow_sum Tot.CC_mean CC.high_mean
## Min. : 0.000 Min. : 0.00000 Min. : 0.00 Min. : 0.00
## 1st Qu.: 0.000 1st Qu.: 0.00000 1st Qu.: 22.20 1st Qu.: 2.20
## Median : 0.100 Median : 0.00000 Median : 52.88 Median :13.23
## Mean : 2.213 Mean : 0.06156 Mean : 50.66 Mean :19.75
## 3rd Qu.: 2.500 3rd Qu.: 0.00000 3rd Qu.: 78.34 3rd Qu.:32.59
## Max. :40.300 Max. :13.86000 Max. :100.00 Max. :94.54
## CC.mid_mean CC.low_mean Sunshine_sum Radiation_sum
## Min. : 0.000 Min. : 0.00 Min. : 0.0 Min. : 283.9
## 1st Qu.: 2.678 1st Qu.: 8.49 1st Qu.: 118.4 1st Qu.:2059.5
## Median : 23.770 Median : 38.70 Median : 371.5 Median :3576.0
## Mean : 31.284 Mean : 39.68 Mean : 377.4 Mean :3945.7
## 3rd Qu.: 54.710 3rd Qu.: 63.58 3rd Qu.: 591.9 3rd Qu.:5730.9
## Max. :100.000 Max. :100.00 Max. :1022.2 Max. :8337.5
## WS.10m_mean WDir.10m_mean WS.80m_mean WDir.80m_mean
## Min. : 2.020 Min. : 34.51 Min. : 2.110 Min. : 36.64
## 1st Qu.: 6.567 1st Qu.:155.87 1st Qu.: 8.738 1st Qu.:158.84
## Median : 9.280 Median :208.99 Median :12.375 Median :215.37
## Mean :10.869 Mean :202.09 Mean :14.393 Mean :207.04
## 3rd Qu.:13.400 3rd Qu.:251.75 3rd Qu.:17.887 3rd Qu.:257.23
## Max. :43.750 Max. :337.65 Max. :54.620 Max. :342.42
## WS.900mb_mean WDir.900mb_mean WGust_mean Temp_max
## Min. : 2.38 Min. : 29.22 Min. : 2.710 Min. :-4.39
## 1st Qu.:13.39 1st Qu.:149.22 1st Qu.: 9.643 1st Qu.:10.65
## Median :20.50 Median :237.40 Median :13.835 Median :16.70
## Mean :25.05 Mean :208.69 Mean :16.791 Mean :16.48
## 3rd Qu.:31.87 3rd Qu.:265.66 3rd Qu.:21.210 3rd Qu.:22.45
## Max. :97.69 Max. :343.97 Max. :70.020 Max. :36.22
## Temp_min Humidity_max Humidity_min MSL_max
## Min. :-11.650 Min. : 55.0 Min. :21.0 Min. : 981.4
## 1st Qu.: 2.998 1st Qu.: 83.0 1st Qu.:45.0 1st Qu.:1015.7
## Median : 8.255 Median : 89.0 Median :54.0 Median :1019.5
## Mean : 8.051 Mean : 87.9 Mean :54.2 Mean :1019.8
## 3rd Qu.: 13.072 3rd Qu.: 94.0 3rd Qu.:63.0 3rd Qu.:1024.4
## Max. : 23.360 Max. :100.0 Max. :93.0 Max. :1041.8
## MSL_min Tot.CC_max Tot.CC_min CC.high_max
## Min. : 976.1 Min. : 0.00 Min. : 0.0 Min. : 0.00
## 1st Qu.:1009.4 1st Qu.:100.00 1st Qu.: 0.0 1st Qu.: 18.00
## Median :1014.4 Median :100.00 Median : 0.0 Median :100.00
## Mean :1014.2 Mean : 89.08 Mean : 8.5 Mean : 61.99
## 3rd Qu.:1019.3 3rd Qu.:100.00 3rd Qu.: 2.1 3rd Qu.:100.00
## Max. :1037.4 Max. :100.00 Max. :100.0 Max. :100.00
## CC.high_min CC.mid_max CC.mid_min CC.low_max
## Min. : 0.0000 Min. : 0.00 Min. : 0.000 Min. : 0.00
## 1st Qu.: 0.0000 1st Qu.: 34.75 1st Qu.: 0.000 1st Qu.:100.00
## Median : 0.0000 Median :100.00 Median : 0.000 Median :100.00
## Mean : 0.7146 Mean : 72.83 Mean : 1.953 Mean : 80.85
## 3rd Qu.: 0.0000 3rd Qu.:100.00 3rd Qu.: 0.000 3rd Qu.:100.00
## Max. :45.0000 Max. :100.00 Max. :100.000 Max. :100.00
## CC.low_min WS.10m_max WS.10m_min WS.80m_max
## Min. : 0.000 Min. : 4.61 Min. : 0.000 Min. : 4.84
## 1st Qu.: 0.000 1st Qu.:12.31 1st Qu.: 1.080 1st Qu.:18.19
## Median : 0.000 Median :17.33 Median : 2.310 Median :23.82
## Mean : 4.551 Mean :19.22 Mean : 3.681 Mean :25.45
## 3rd Qu.: 0.000 3rd Qu.:23.88 3rd Qu.: 4.800 3rd Qu.:29.92
## Max. :100.000 Max. :69.42 Max. :36.040 Max. :79.78
## WS.80m_min WS.900mb_max WS.900mb_min WGust_max
## Min. : 0.000 Min. : 6.13 Min. : 0.00 Min. : 5.04
## 1st Qu.: 1.125 1st Qu.: 25.18 1st Qu.: 3.26 1st Qu.:19.08
## Median : 2.550 Median : 37.94 Median : 7.15 Median :25.92
## Mean : 4.901 Mean : 42.26 Mean :11.37 Mean :29.49
## 3rd Qu.: 6.130 3rd Qu.: 55.51 3rd Qu.:14.69 3rd Qu.:36.36
## Max. :46.130 Max. :124.20 Max. :82.07 Max. :97.92
## WGust_min pluie.demain
## Min. : 0.000 Mode :logical
## 1st Qu.: 2.160 FALSE:625
## Median : 3.960 TRUE :619
## Mean : 6.655
## 3rd Qu.: 8.640
## Max. :51.120On remarque qu’il pleut 1 jour sur deux en moyenne.
Les variable Hour et Minutes valent toutes 0 et peuvent donc être éliminées directement.
2. Construction du modèle saturé
On est dans le cadre d’un modèle logistique simple. La variable à estimer, pluie.demain est une variable binaire 0-1. Et suit une loi de Bernoulli de probabilité approximative 1/2. Pour estimer la variable pluie.demain on va donc effectuer une regression logistique sous R en utilisant la fonction glm(). Par défaut glm utilise comme fonction de lien la fonction logit. Les deux modèles extrêmes sont les modèles saturé et le modèle null (où toutes les coefficients des covraiables sont nulles). On ne s’intéresse ici qu’au modèle saturé qui va nous permetrre de construire nos modèles intermédiares à k<44 covariables.
logit_sat = glm(pluie.demain ~ . , family = binomial, data = meteoTrain, singular.ok = TRUE )Estimateur des parametres du modèle -> coef
Aucuns des coéfficients n’est à NA, toutes les variables sont bien interprétées.
coef(logit_sat)[is.na(coef(logit_sat))]## named numeric(0)table(meteoTrain$Hour,meteoTrain$Minute)## < table of extent 0 x 0 >Anova
anova(logit_sat,test="LR")## Analysis of Deviance Table
##
## Model: binomial, link: logit
##
## Response: pluie.demain
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 1243 1724.5
## X 1 0.021 1242 1724.5 0.8857467
## Year 1 1.717 1241 1722.8 0.1900469
## Month 1 1.411 1240 1721.4 0.2349204
## Day 1 0.006 1239 1721.4 0.9370321
## Temp_mean 1 3.795 1238 1717.6 0.0513973 .
## Humidity_mean 1 10.390 1237 1707.2 0.0012673 **
## MSL_mean 1 185.428 1236 1521.8 < 2.2e-16 ***
## Tot.Prcp_sum 1 13.074 1235 1508.7 0.0002993 ***
## Snow_sum 1 3.214 1234 1505.5 0.0730026 .
## Tot.CC_mean 1 27.870 1233 1477.6 1.298e-07 ***
## CC.high_mean 1 14.417 1232 1463.2 0.0001465 ***
## CC.mid_mean 1 2.729 1231 1460.5 0.0985610 .
## CC.low_mean 1 0.003 1230 1460.5 0.9558099
## Sunshine_sum 1 0.049 1229 1460.4 0.8248690
## Radiation_sum 1 0.007 1228 1460.4 0.9314771
## WS.10m_mean 1 0.227 1227 1460.2 0.6336471
## WDir.10m_mean 1 2.030 1226 1458.1 0.1542416
## WS.80m_mean 1 1.911 1225 1456.2 0.1668453
## WDir.80m_mean 1 1.138 1224 1455.1 0.2861729
## WS.900mb_mean 1 3.207 1223 1451.9 0.0733073 .
## WDir.900mb_mean 1 28.340 1222 1423.5 1.018e-07 ***
## WGust_mean 1 1.890 1221 1421.7 0.1691643
## Temp_max 1 6.639 1220 1415.0 0.0099768 **
## Temp_min 1 8.776 1219 1406.2 0.0030514 **
## Humidity_max 1 0.056 1218 1406.2 0.8122857
## Humidity_min 1 1.543 1217 1404.6 0.2142302
## MSL_max 1 0.994 1216 1403.6 0.3188399
## MSL_min 1 6.260 1215 1397.4 0.0123496 *
## Tot.CC_max 1 0.539 1214 1396.8 0.4629026
## Tot.CC_min 1 0.001 1213 1396.8 0.9767642
## CC.high_max 1 7.464 1212 1389.4 0.0062939 **
## CC.high_min 1 0.931 1211 1388.4 0.3346976
## CC.mid_max 1 12.976 1210 1375.5 0.0003155 ***
## CC.mid_min 1 1.540 1209 1373.9 0.2146728
## CC.low_max 1 0.004 1208 1373.9 0.9494757
## CC.low_min 1 5.470 1207 1368.5 0.0193412 *
## WS.10m_max 1 0.895 1206 1367.6 0.3441553
## WS.10m_min 1 7.579 1205 1360.0 0.0059044 **
## WS.80m_max 1 0.364 1204 1359.6 0.5461963
## WS.80m_min 1 0.004 1203 1359.6 0.9497787
## WS.900mb_max 1 0.866 1202 1358.8 0.3520559
## WS.900mb_min 1 1.586 1201 1357.2 0.2079170
## WGust_max 1 0.792 1200 1356.4 0.3736023
## WGust_min 1 0.498 1199 1355.9 0.4803746
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1anova(logit_sat,test="Cp")## Analysis of Deviance Table
##
## Model: binomial, link: logit
##
## Response: pluie.demain
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Cp
## NULL 1243 1724.5 1726.5
## X 1 0.021 1242 1724.5 1728.5
## Year 1 1.717 1241 1722.8 1728.8
## Month 1 1.411 1240 1721.4 1729.4
## Day 1 0.006 1239 1721.4 1731.4
## Temp_mean 1 3.795 1238 1717.6 1729.6
## Humidity_mean 1 10.390 1237 1707.2 1721.2
## MSL_mean 1 185.428 1236 1521.8 1537.8
## Tot.Prcp_sum 1 13.074 1235 1508.7 1526.7
## Snow_sum 1 3.214 1234 1505.5 1525.5
## Tot.CC_mean 1 27.870 1233 1477.6 1499.6
## CC.high_mean 1 14.417 1232 1463.2 1487.2
## CC.mid_mean 1 2.729 1231 1460.5 1486.4
## CC.low_mean 1 0.003 1230 1460.5 1488.4
## Sunshine_sum 1 0.049 1229 1460.4 1490.4
## Radiation_sum 1 0.007 1228 1460.4 1492.4
## WS.10m_mean 1 0.227 1227 1460.2 1494.2
## WDir.10m_mean 1 2.030 1226 1458.1 1494.1
## WS.80m_mean 1 1.911 1225 1456.2 1494.2
## WDir.80m_mean 1 1.138 1224 1455.1 1495.1
## WS.900mb_mean 1 3.207 1223 1451.9 1493.9
## WDir.900mb_mean 1 28.340 1222 1423.5 1467.5
## WGust_mean 1 1.890 1221 1421.7 1467.6
## Temp_max 1 6.639 1220 1415.0 1463.0
## Temp_min 1 8.776 1219 1406.2 1456.2
## Humidity_max 1 0.056 1218 1406.2 1458.2
## Humidity_min 1 1.543 1217 1404.6 1458.6
## MSL_max 1 0.994 1216 1403.6 1459.6
## MSL_min 1 6.260 1215 1397.4 1455.4
## Tot.CC_max 1 0.539 1214 1396.8 1456.8
## Tot.CC_min 1 0.001 1213 1396.8 1458.8
## CC.high_max 1 7.464 1212 1389.4 1453.4
## CC.high_min 1 0.931 1211 1388.4 1454.4
## CC.mid_max 1 12.976 1210 1375.5 1443.5
## CC.mid_min 1 1.540 1209 1373.9 1443.9
## CC.low_max 1 0.004 1208 1373.9 1445.9
## CC.low_min 1 5.470 1207 1368.5 1442.5
## WS.10m_max 1 0.895 1206 1367.6 1443.6
## WS.10m_min 1 7.579 1205 1360.0 1438.0
## WS.80m_max 1 0.364 1204 1359.6 1439.6
## WS.80m_min 1 0.004 1203 1359.6 1441.6
## WS.900mb_max 1 0.866 1202 1358.8 1442.7
## WS.900mb_min 1 1.586 1201 1357.2 1443.2
## WGust_max 1 0.792 1200 1356.4 1444.4
## WGust_min 1 0.498 1199 1355.9 1445.9On obtient le premier modèle “saturé”: m_sature
m_sature = glm(formula = pluie.demain ~ (X+Year+Month+Day+
Temp_mean + Humidity_mean + MSL_mean + Tot.Prcp_sum + Snow_sum +
Tot.CC_mean + CC.high_mean + CC.mid_mean + CC.low_mean +
Sunshine_sum + Radiation_sum + WS.10m_mean + WDir.10m_mean +
WS.80m_mean + WDir.80m_mean + WS.900mb_mean + WDir.900mb_mean +
WGust_mean + Temp_max + Temp_min + Humidity_max + Humidity_min +
MSL_max + MSL_min + Tot.CC_max + Tot.CC_min + CC.high_max +
CC.high_min + CC.mid_max + CC.mid_min + CC.low_max + CC.low_min +
WS.10m_max + WS.10m_min + WS.80m_max + WS.80m_min + WS.900mb_max +
WS.900mb_min + WGust_max + WGust_min), family = binomial, data = meteoTrain)
summary(m_sature)##
## Call:
## glm(formula = pluie.demain ~ (X + Year + Month + Day + Temp_mean +
## Humidity_mean + MSL_mean + Tot.Prcp_sum + Snow_sum + Tot.CC_mean +
## CC.high_mean + CC.mid_mean + CC.low_mean + Sunshine_sum +
## Radiation_sum + WS.10m_mean + WDir.10m_mean + WS.80m_mean +
## WDir.80m_mean + WS.900mb_mean + WDir.900mb_mean + WGust_mean +
## Temp_max + Temp_min + Humidity_max + Humidity_min + MSL_max +
## MSL_min + Tot.CC_max + Tot.CC_min + CC.high_max + CC.high_min +
## CC.mid_max + CC.mid_min + CC.low_max + CC.low_min + WS.10m_max +
## WS.10m_min + WS.80m_max + WS.80m_min + WS.900mb_max + WS.900mb_min +
## WGust_max + WGust_min), family = binomial, data = meteoTrain)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.7555 -0.8855 -0.2028 0.8975 2.5328
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.157e+04 8.898e+04 -0.130 0.896503
## X -1.576e-02 1.212e-01 -0.130 0.896512
## Year 5.800e+00 4.426e+01 0.131 0.895741
## Month 4.784e-01 3.686e+00 0.130 0.896716
## Day 7.743e-03 1.211e-01 0.064 0.949024
## Temp_mean 1.423e-01 1.613e-01 0.882 0.377716
## Humidity_mean -3.554e-02 2.982e-02 -1.192 0.233353
## MSL_mean 1.558e-01 1.249e-01 1.247 0.212493
## Tot.Prcp_sum 1.791e-02 2.408e-02 0.744 0.457175
## Snow_sum 4.245e-01 3.626e-01 1.171 0.241691
## Tot.CC_mean -1.842e-03 1.199e-02 -0.154 0.877932
## CC.high_mean 4.772e-03 6.748e-03 0.707 0.479480
## CC.mid_mean -3.438e-03 6.390e-03 -0.538 0.590582
## CC.low_mean 7.678e-03 8.201e-03 0.936 0.349124
## Sunshine_sum -5.701e-04 8.643e-04 -0.660 0.509474
## Radiation_sum 4.491e-05 9.752e-05 0.461 0.645157
## WS.10m_mean -2.544e-02 9.423e-02 -0.270 0.787139
## WDir.10m_mean 7.023e-03 5.080e-03 1.383 0.166794
## WS.80m_mean -1.226e-01 6.524e-02 -1.879 0.060258 .
## WDir.80m_mean -8.996e-03 5.265e-03 -1.709 0.087504 .
## WS.900mb_mean 1.940e-02 2.259e-02 0.859 0.390316
## WDir.900mb_mean 6.983e-03 1.370e-03 5.099 3.42e-07 ***
## WGust_mean 3.250e-02 3.465e-02 0.938 0.348332
## Temp_max 4.438e-02 9.357e-02 0.474 0.635298
## Temp_min -2.069e-01 8.494e-02 -2.436 0.014831 *
## Humidity_max 2.259e-02 1.934e-02 1.168 0.242842
## Humidity_min 2.590e-02 1.638e-02 1.581 0.113914
## MSL_max -1.428e-01 6.924e-02 -2.062 0.039194 *
## MSL_min -1.007e-01 6.541e-02 -1.540 0.123533
## Tot.CC_max -5.037e-03 4.718e-03 -1.068 0.285715
## Tot.CC_min 1.420e-02 7.796e-03 1.821 0.068557 .
## CC.high_max 2.525e-03 2.738e-03 0.922 0.356473
## CC.high_min -2.258e-02 2.568e-02 -0.880 0.379088
## CC.mid_max 1.085e-02 3.070e-03 3.534 0.000409 ***
## CC.mid_min 6.689e-03 9.757e-03 0.686 0.493004
## CC.low_max -5.397e-04 3.321e-03 -0.163 0.870873
## CC.low_min -1.640e-02 7.881e-03 -2.081 0.037405 *
## WS.10m_max 2.694e-02 3.303e-02 0.816 0.414618
## WS.10m_min 1.095e-01 5.766e-02 1.900 0.057440 .
## WS.80m_max 2.306e-03 2.711e-02 0.085 0.932223
## WS.80m_min 1.091e-02 4.028e-02 0.271 0.786559
## WS.900mb_max 1.473e-03 1.070e-02 0.138 0.890472
## WS.900mb_min -1.657e-02 1.696e-02 -0.977 0.328711
## WGust_max 1.037e-02 1.693e-02 0.612 0.540299
## WGust_min -1.831e-02 2.596e-02 -0.705 0.480510
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1724.5 on 1243 degrees of freedom
## Residual deviance: 1355.9 on 1199 degrees of freedom
## AIC: 1445.9
##
## Number of Fisher Scoring iterations: 6Points influents / abérrants
which(influence(m_sature)$hat>0.2)## 95 270 716
## 95 270 716## X Year Month Day Temp_mean Humidity_mean MSL_mean Tot.Prcp_sum
## 95 189 2010 12 6 2.61 93.88 999.14 40.3
## Snow_sum Tot.CC_mean CC.high_mean CC.mid_mean CC.low_mean Sunshine_sum
## 95 0 100 91.96 100 100 0
## Radiation_sum WS.10m_mean WDir.10m_mean WS.80m_mean WDir.80m_mean
## 95 283.91 15.89 101.5 20.81 108.91
## WS.900mb_mean WDir.900mb_mean WGust_mean Temp_max Temp_min Humidity_max
## 95 17.1 219.84 14.55 3.59 1.74 98
## Humidity_min MSL_max MSL_min Tot.CC_max Tot.CC_min CC.high_max
## 95 92 1004 996.4 100 100 100
## CC.high_min CC.mid_max CC.mid_min CC.low_max CC.low_min WS.10m_max
## 95 45 100 100 100 100 24.07
## WS.10m_min WS.80m_max WS.80m_min WS.900mb_max WS.900mb_min WGust_max
## 95 3.32 32.21 7.21 33.59 7.99 24.12
## WGust_min pluie.demain
## 95 0.72 TRUE## X Year Month Day Temp_mean Humidity_mean MSL_mean Tot.Prcp_sum
## 270 607 2012 1 28 4.53 79.21 1025.41 0
## Snow_sum Tot.CC_mean CC.high_mean CC.mid_mean CC.low_mean Sunshine_sum
## 270 0 100 20.92 100 95.62 0
## Radiation_sum WS.10m_mean WDir.10m_mean WS.80m_mean WDir.80m_mean
## 270 1346.57 6.62 118.29 8.59 136.23
## WS.900mb_mean WDir.900mb_mean WGust_mean Temp_max Temp_min
## 270 13.53 107.95 10.06 6.09 2.39
## Humidity_max Humidity_min MSL_max MSL_min Tot.CC_max Tot.CC_min
## 270 90 67 1026.5 1024 100 100
## CC.high_max CC.high_min CC.mid_max CC.mid_min CC.low_max CC.low_min
## 270 100 0 100 100 100 7
## WS.10m_max WS.10m_min WS.80m_max WS.80m_min WS.900mb_max WS.900mb_min
## 270 13.39 0.72 19.38 0.36 29.83 2.19
## WGust_max WGust_min pluie.demain
## 270 23.76 0.72 FALSE## X Year Month Day Temp_mean Humidity_mean MSL_mean Tot.Prcp_sum
## 716 1671 2014 12 27 2.3 79.62 1005.77 6.3
## Snow_sum Tot.CC_mean CC.high_mean CC.mid_mean CC.low_mean Sunshine_sum
## 716 1.82 93.42 40.38 86.38 89.83 0
## Radiation_sum WS.10m_mean WDir.10m_mean WS.80m_mean WDir.80m_mean
## 716 1092.03 42.31 213.72 50.88 216.3
## WS.900mb_mean WDir.900mb_mean WGust_mean Temp_max Temp_min
## 716 75.5 252.54 48.04 5.33 -0.52
## Humidity_max Humidity_min MSL_max MSL_min Tot.CC_max Tot.CC_min
## 716 89 62 1025.7 996 100 7
## CC.high_max CC.high_min CC.mid_max CC.mid_min CC.low_max CC.low_min
## 716 100 0 100 0 100 0
## WS.10m_max WS.10m_min WS.80m_max WS.80m_min WS.900mb_max WS.900mb_min
## 716 69.42 6.62 79.78 9.18 109.18 16.3
## WGust_max WGust_min pluie.demain
## 716 77.4 3.6 TRUE
which(cooks.distance(m_sature)>0.03)## 397
## 397Il conviendrait d’analyser pourquoi ces jours sont si particuliers et l’élimintaion direste n’est pas conseillée. Cependant vu l’e nombre de jours’échantillon relativement important, on va tout de même les supprimer. Sauvegarde des données et suppression des lignes 95,270,716 et 397
Rechargement du modèle saturé:
m_sature = glm(formula = pluie.demain ~ (X+Year+Month+Day+
Temp_mean + Humidity_mean + MSL_mean + Tot.Prcp_sum + Snow_sum +
Tot.CC_mean + CC.high_mean + CC.mid_mean + CC.low_mean +
Sunshine_sum + Radiation_sum + WS.10m_mean + WDir.10m_mean +
WS.80m_mean + WDir.80m_mean + WS.900mb_mean + WDir.900mb_mean +
WGust_mean + Temp_max + Temp_min + Humidity_max + Humidity_min +
MSL_max + MSL_min + Tot.CC_max + Tot.CC_min + CC.high_max +
CC.high_min + CC.mid_max + CC.mid_min + CC.low_max + CC.low_min +
WS.10m_max + WS.10m_min + WS.80m_max + WS.80m_min + WS.900mb_max +
WS.900mb_min + WGust_max + WGust_min), family = binomial, data = meteoTrain)Détection de colinéarité
Les corrélations positives sont affichées en bleu et les corrélations négatives en rouge. L’intensité de la couleur et la taille des cercles sont proportionnelles aux coefficients de corrélation.
On retrouve des corrélations fortes entre le groupe de vriables de type vents/vitesse/rafale. Ce qui est normale et attendu. Entre nébulosité et rayonnement une corrélation négative importante. Mais dans l’ensemble et entre nos différents types de variables, il y a finalement assez peu de corélations. Ce qui est une bonne chose. Notament pour utiliser des tests anova entre modèles emboîtés.
Etude des variables temporelles: Day, Month, Year, X
logit_sat_Day = glm(pluie.demain ~ . -Day , family = binomial, data = meteoTrain )
logit_sat_DayMonth = glm(pluie.demain ~ . -Day -Month, family = binomial, data = meteoTrain )
logit_sat_DayMonthYear = glm(pluie.demain ~ . -Day -Month -Year, family = binomial, data = meteoTrain )
logit_sat_DayMonthYearX = glm(pluie.demain ~ . -Day -Month -Year -X, family = binomial, data = meteoTrain )Annova, test de déviance entre modèles emboîtés
anova(m_sature,logit_sat_Day,test="LRT")
anova(m_sature,logit_sat_DayMonth,test="LRT")
anova(m_sature,logit_sat_DayMonthYear,test="LRT")
anova(m_sature,logit_sat_DayMonthYearX,test="LRT")Les tests anova successifs sur les modèles emboîtés, où l’on supprime à chaque étape l’une des variables temporelles: Day, Month, Year, X montrent que l’on peut se passer de ces variables. la p-value (Pr(>Chi)) est importante > 0.6 dans tous les cas. On conserve donc l’hypothèse H0 de nullité des coéfficients, c’est à dire le plus petit modèle.
On peut étudier et faire une première analyse des données : PCA en utilsant le package FactoMineR. Et essayer à partir de cette analyse de sélectionner des variables en fonction de leur bonne contribution et représentation aux axes principaux.
3. Analyse PCA
Les deux variable Year et X (ordonnancement) sont passée en variables quanitatives supplémentaires. Ce ne sont pas des variables météorologiques. Après une première étude, ces variables perturbent plus l’analyse qu’elles n’apportent d’information. En particulier au niveau de l’analyse du 5ème axe. On a vu aussi en fin du chapitre précédent §2, avec des test anova, que les variables temporelles ne semblent pas sigificatives. Ce ne sont pas aussi des variables météorologiques, on va aussi les passer en variables quantitatives supplémentaires.
library(readxl)
library(VIM)
library(FactoMineR)
library(factoextra)
library(corrplot)
res.pca.meteoTrain<-PCA(meteoTrain,graph=FALSE, quali.sup=45, quanti.sup = 1:4)
summary(res.pca.meteoTrain)##
## Call:
## PCA(X = meteoTrain, quanti.sup = 1:4, quali.sup = 45, graph = FALSE)
##
##
## Eigenvalues
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5 Dim.6
## Variance 12.994 5.449 4.382 2.404 1.967 1.897
## % of var. 32.485 13.623 10.955 6.010 4.918 4.744
## Cumulative % of var. 32.485 46.108 57.063 63.073 67.991 72.735
## Dim.7 Dim.8 Dim.9 Dim.10 Dim.11 Dim.12
## Variance 1.332 1.198 1.116 0.941 0.920 0.683
## % of var. 3.330 2.996 2.790 2.353 2.299 1.706
## Cumulative % of var. 76.065 79.061 81.851 84.203 86.503 88.209
## Dim.13 Dim.14 Dim.15 Dim.16 Dim.17 Dim.18
## Variance 0.633 0.577 0.516 0.471 0.371 0.351
## % of var. 1.582 1.441 1.289 1.177 0.929 0.877
## Cumulative % of var. 89.791 91.232 92.521 93.699 94.627 95.504
## Dim.19 Dim.20 Dim.21 Dim.22 Dim.23 Dim.24
## Variance 0.228 0.201 0.187 0.158 0.146 0.141
## % of var. 0.569 0.502 0.468 0.395 0.365 0.351
## Cumulative % of var. 96.073 96.575 97.044 97.439 97.804 98.155
## Dim.25 Dim.26 Dim.27 Dim.28 Dim.29 Dim.30
## Variance 0.117 0.112 0.102 0.068 0.061 0.058
## % of var. 0.293 0.281 0.255 0.171 0.153 0.144
## Cumulative % of var. 98.448 98.729 98.985 99.155 99.308 99.452
## Dim.31 Dim.32 Dim.33 Dim.34 Dim.35 Dim.36
## Variance 0.047 0.041 0.031 0.025 0.023 0.020
## % of var. 0.117 0.102 0.077 0.063 0.057 0.050
## Cumulative % of var. 99.569 99.671 99.748 99.812 99.869 99.918
## Dim.37 Dim.38 Dim.39 Dim.40
## Variance 0.019 0.008 0.003 0.002
## % of var. 0.048 0.020 0.008 0.006
## Cumulative % of var. 99.966 99.986 99.994 100.000
##
## Individuals (the 10 first)
## Dist Dim.1 ctr cos2 Dim.2 ctr cos2
## 1 | 3.691 | -0.297 0.001 0.006 | 2.319 0.080 0.395 |
## 2 | 4.155 | 0.100 0.000 0.001 | 2.444 0.088 0.346 |
## 3 | 6.939 | -5.804 0.209 0.700 | -1.643 0.040 0.056 |
## 4 | 3.934 | -0.910 0.005 0.053 | 1.321 0.026 0.113 |
## 5 | 5.052 | -0.324 0.001 0.004 | 1.660 0.041 0.108 |
## 6 | 4.835 | -1.025 0.007 0.045 | 1.171 0.020 0.059 |
## 7 | 5.451 | 1.099 0.007 0.041 | 3.901 0.225 0.512 |
## 8 | 4.996 | 0.950 0.006 0.036 | 2.919 0.126 0.341 |
## 9 | 8.644 | 2.257 0.032 0.068 | 5.855 0.507 0.459 |
## 10 | 5.308 | 2.801 0.049 0.278 | 2.653 0.104 0.250 |
## Dim.3 ctr cos2
## 1 0.431 0.003 0.014 |
## 2 -0.719 0.010 0.030 |
## 3 0.717 0.009 0.011 |
## 4 2.231 0.092 0.322 |
## 5 3.694 0.251 0.535 |
## 6 1.694 0.053 0.123 |
## 7 -0.118 0.000 0.000 |
## 8 -0.063 0.000 0.000 |
## 9 -0.827 0.013 0.009 |
## 10 1.025 0.019 0.037 |
##
## Variables (the 10 first)
## Dim.1 ctr cos2 Dim.2 ctr cos2 Dim.3
## Temp_mean | -0.398 1.219 0.158 | 0.070 0.091 0.005 | 0.778
## Humidity_mean | 0.411 1.301 0.169 | 0.494 4.476 0.244 | -0.521
## MSL_mean | -0.410 1.294 0.168 | -0.217 0.863 0.047 | -0.592
## Tot.Prcp_sum | 0.489 1.839 0.239 | 0.340 2.126 0.116 | 0.150
## Snow_sum | 0.162 0.202 0.026 | 0.074 0.101 0.005 | -0.101
## Tot.CC_mean | 0.741 4.224 0.549 | 0.542 5.401 0.294 | -0.008
## CC.high_mean | 0.384 1.135 0.148 | 0.431 3.412 0.186 | 0.339
## CC.mid_mean | 0.655 3.303 0.429 | 0.476 4.156 0.226 | 0.252
## CC.low_mean | 0.693 3.695 0.480 | 0.470 4.056 0.221 | -0.152
## Sunshine_sum | -0.727 4.066 0.528 | -0.468 4.015 0.219 | 0.191
## ctr cos2
## Temp_mean 13.815 0.605 |
## Humidity_mean 6.193 0.271 |
## MSL_mean 8.002 0.351 |
## Tot.Prcp_sum 0.511 0.022 |
## Snow_sum 0.233 0.010 |
## Tot.CC_mean 0.001 0.000 |
## CC.high_mean 2.618 0.115 |
## CC.mid_mean 1.454 0.064 |
## CC.low_mean 0.528 0.023 |
## Sunshine_sum 0.828 0.036 |
##
## Supplementary continuous variables
## Dim.1 cos2 Dim.2 cos2 Dim.3 cos2
## X | 0.026 0.001 | -0.117 0.014 | 0.062 0.004 |
## Year | 0.031 0.001 | -0.122 0.015 | 0.056 0.003 |
## Month | -0.049 0.002 | 0.050 0.003 | 0.040 0.002 |
## Day | -0.016 0.000 | 0.022 0.000 | -0.012 0.000 |
##
## Supplementary categories
## Dist Dim.1 cos2 v.test Dim.2 cos2
## FALSE | 1.282 | -1.055 0.678 -10.351 | -0.407 0.101
## TRUE | 1.294 | 1.065 0.678 10.351 | 0.411 0.101
## v.test Dim.3 cos2 v.test
## FALSE -6.165 | -0.558 0.189 -9.427 |
## TRUE 6.165 | 0.563 0.189 9.427 |3.1 Choix du nombre d’axes
On cherche une cassure/décroissance au niveau des inerties. Une première cassure dès la 1ere dimension, puis après les 3 premières dimensions. puis une 3ème après la 6ème dimension. Et une 4ème au nivau de la 11ème dimension. On entrevoit des difficultés à réduire le nombre de dimensions. La décroissance est lente.
Il va être difficile de restreindre le nombre des variables explicatives dans notre modèle glm.
Décomposition de la variabilité par axes:
Pourcentage de l’inertie par rapport aux dimensions (axes).Qui mesure la qualité de représentation des données / axes. Ici le 1er axe représente 32% de la variabilité des données. Et c’est aussi une mesure de l’importance relative des axes.Ici le 1er exprime 2.5 fois plus que le 2ème.
La décroissance est forte pour les 3 -1ères valeur propres. Puis de l’ordre de 5% pour les 4 suivantes. Et ensuite la décroissance est régulière.
- les 3 -1ères dimensions expliquent 57% de la variabilité.
- les 5 -1ères dimensions expliquent 68% de la variabilité.
- les 9 -1ères dimensions expliquent 80% de la variabilité.
- les 14 -1ères dimensions expliquent 91% de la variabilité.
Critère du choix du nombre d’axes:
Critère du coude atteint au 1er ou 3éme axe. Pour le critère de kaiser on regarde la 2ème composante et on compare : variance cumulative >= I/p = 100/40 = 2.5 soit à partir de comp 9. Pour le critère des 80% on regarde la 3ème composante: atteint à partir de la dimension 9.
On remarque qu’il est difficle de réduire la dimension. On va se contenter d’éudier les 5 premiers axes qui expliquent 68% de l’inertie.
## eigenvalue percentage of variance
## comp 1 12.99 32.49
## comp 2 5.45 13.62
## comp 3 4.38 10.95
## comp 4 2.40 6.01
## comp 5 1.97 4.92
## comp 6 1.90 4.74
## comp 7 1.33 3.33
## comp 8 1.20 3.00
## comp 9 1.12 2.79
## comp 10 0.94 2.35
## comp 11 0.92 2.30
## comp 12 0.68 1.71
## comp 13 0.63 1.58
## comp 14 0.58 1.44
## comp 15 0.52 1.29
## comp 16 0.47 1.18
## cumulative percentage of variance
## comp 1 32.49
## comp 2 46.11
## comp 3 57.06
## comp 4 63.07
## comp 5 67.99
## comp 6 72.73
## comp 7 76.06
## comp 8 79.06
## comp 9 81.85
## comp 10 84.20
## comp 11 86.50
## comp 12 88.21
## comp 13 89.79
## comp 14 91.23
## comp 15 92.52
## comp 16 93.703.2 Analyse des résultats - Graphique des individus et des variables
3.2.1 Projection du nuage des individus (jours de pluie/sec) sur les 2 axes principaux
Le plan principale représente 46 % environ de l’inertie totale. Le nuage a une allure régulière. Les données Pluie/sec sont relativement centrées / axe1 et axe2.
3.2.2 Projection des variables sur les 2 axes principaux
Les jours secs ont leur éllipse de confiance décalé dans la direction du vecteur “Sunshine Sum”
L’axe 2 oppose les vriables Nébulosité (Cloud Cover) à Sunshine Sum (Ensoleillement)
Quart Bas/Droit: Contributions > 0 avec la 1ere dimension. On retrouve en majorité des données de type Wind Speed et Wind Gust Données très corrélées entre elles et avec le 1er axes. Ceux sont ces variables (de type Wind-Speed/Wind-Gust) qui contribue le plus à l’axe 1.
Quart Haut/Droit:On trouve en majorité des données de catégories Cloud Cover et Wind direction (de manière moins évidentes) Aussi Humidity. Contributions principales >0 au 2ème axe.
Quart Haut/Gauche: Données de type température. Contribution faible à l’une ou l’autre des dimension
Quart Bas/Gauche: Données radiation et sunshine duration en corrélation négative En opposition avec les données de type Cloud Coverage. En effet plus la couverture nuageuse est importante et plus courte est la durée du coucher de soleil ou bien la quantité de radiation.
3.2.3 Autres axes
Projection du nuage des individus (journées) sur les axes 2 et 3
- Le nuage a une allure régulière. Les données Pluie/sec sont relativement centrées / axe2 et axe3. Etalement et applatissement le long de l’axe 3.
Les variables de type température sont parfaitement corrélées avec l’axe 3. En opposition avec la Pression Max: MSL_Max
Projection du nuage des individus (journées) sur les axes 3 et 4
- Le nuage a une allure régulière. Les données Pluie/sec sont relativement centrées / axe3 et axe4.
- Pas de remarques paticulières. Les nuages sont toujours centrés. Le décalage pluie/sec est faible.
On va maintenant passer aux analyses automatiques et détaillées de chacun des axes. Et étudier la contribution et représentation des variables.
3.3 Analyses automatiques
3.3.1 coefficients de corrélation entre chacunes des variables et les 5 premières composantes principales
(ce qui correspond aux coordonnées des individus sur les 5 premiers axes)
A partir de la fonction dimdesc
dimdesc <- dimdesc(res.pca.meteoTrain, proba=1e-5) #, nbelements = 16)
dimdesc## $Dim.1
## $Dim.1$quanti
## correlation p.value
## WGust_mean 0.8574807 0.000000e+00
## WS.10m_mean 0.8440876 0.000000e+00
## WS.80m_mean 0.8356130 0.000000e+00
## WS.10m_max 0.8165994 9.327448e-298
## WGust_max 0.8099670 3.468518e-289
## WS.80m_max 0.8000925 4.930158e-277
## WS.900mb_mean 0.7963525 1.314316e-272
## WS.900mb_max 0.7719127 8.280547e-246
## Tot.CC_mean 0.7408979 2.901487e-216
## WGust_min 0.7136223 1.382084e-193
## CC.low_mean 0.6928993 4.612783e-178
## WS.900mb_min 0.6913563 5.849040e-177
## WS.80m_min 0.6853018 1.069607e-172
## WS.10m_min 0.6754141 5.891137e-166
## CC.mid_mean 0.6550850 6.796469e-153
## CC.low_max 0.5217305 1.665658e-87
## CC.mid_max 0.5039330 7.870246e-81
## Tot.Prcp_sum 0.4889016 1.712046e-75
## Tot.CC_max 0.4578505 2.828133e-65
## Humidity_min 0.4412918 3.049997e-60
## Humidity_mean 0.4110862 9.452170e-52
## Tot.CC_min 0.3864767 1.895459e-45
## WDir.80m_mean 0.3857604 2.839632e-45
## WDir.10m_mean 0.3856171 3.078358e-45
## CC.high_mean 0.3841131 7.166886e-45
## WDir.900mb_mean 0.3664061 1.081889e-40
## CC.high_max 0.3232119 1.520116e-31
## CC.low_min 0.2749393 6.050828e-23
## CC.mid_min 0.2665186 1.311338e-21
## CC.high_min 0.1798096 1.801093e-10
## Humidity_max 0.1787164 2.327696e-10
## Snow_sum 0.1620021 9.644545e-09
## Temp_min -0.2866972 6.865378e-25
## MSL_max -0.3348556 7.189327e-34
## Temp_mean -0.3980382 2.416504e-48
## MSL_mean -0.4101024 1.728873e-51
## MSL_min -0.4607706 3.426275e-66
## Temp_max -0.4747598 1.042050e-70
## Radiation_sum -0.6228247 3.971497e-134
## Sunshine_sum -0.7268956 2.811213e-204
##
## $Dim.1$quali
## R2 p.value
## pluie.demain 0.08648185 3.692656e-26
##
## $Dim.1$category
## Estimate p.value
## pluie.demain=TRUE 1.060082 3.692656e-26
## pluie.demain=FALSE -1.060082 3.692656e-26
##
##
## $Dim.2
## $Dim.2$quanti
## correlation p.value
## Tot.CC_mean 0.5424882 8.290826e-96
## Humidity_mean 0.4938791 3.129205e-77
## CC.mid_mean 0.4758820 4.429442e-71
## CC.low_mean 0.4701351 3.423062e-69
## Tot.CC_min 0.4541764 3.910910e-64
## Humidity_min 0.4486877 1.865118e-62
## CC.high_mean 0.4311984 2.619641e-57
## CC.high_max 0.4228479 5.908875e-55
## Humidity_max 0.4204956 2.645563e-54
## Tot.CC_max 0.4098892 1.970088e-51
## CC.mid_max 0.3953804 1.145039e-47
## CC.low_max 0.3805095 5.332482e-44
## CC.low_min 0.3571271 1.323932e-38
## Tot.Prcp_sum 0.3404032 5.164506e-35
## WDir.80m_mean 0.2810643 6.028830e-24
## WDir.10m_mean 0.2630468 4.516601e-21
## CC.mid_min 0.2529255 1.497977e-19
## Temp_min 0.1523154 7.059711e-08
## CC.high_min 0.1387599 9.339973e-07
## WDir.900mb_mean 0.1358303 1.582551e-06
## MSL_min -0.1838703 6.849745e-11
## Radiation_sum -0.2072085 1.714905e-13
## MSL_mean -0.2168806 1.149361e-14
## MSL_max -0.2366163 3.069114e-17
## WGust_max -0.3115813 2.539498e-29
## WS.900mb_max -0.3216959 3.000732e-31
## WS.10m_max -0.3707753 1.064362e-41
## WS.80m_max -0.3796363 8.639317e-44
## WS.900mb_mean -0.4093362 2.763117e-51
## WGust_mean -0.4189961 6.836640e-54
## WS.10m_mean -0.4526597 1.145925e-63
## WGust_min -0.4602742 4.912432e-66
## WS.900mb_min -0.4629140 7.182883e-67
## WS.80m_mean -0.4643733 2.463491e-67
## Sunshine_sum -0.4677459 2.036156e-68
## WS.10m_min -0.4728323 4.495346e-70
## WS.80m_min -0.4733857 2.957406e-70
##
## $Dim.2$quali
## R2 p.value
## pluie.demain 0.03067905 5.309807e-10
##
## $Dim.2$category
## Estimate p.value
## pluie.demain=TRUE 0.4088742 5.309807e-10
## pluie.demain=FALSE -0.4088742 5.309807e-10
##
##
## $Dim.3
## $Dim.3$quanti
## correlation p.value
## Temp_mean 0.7780584 3.181012e-252
## Temp_max 0.7635228 2.279685e-237
## Temp_min 0.7513588 9.735593e-226
## Radiation_sum 0.5338218 2.854514e-92
## CC.high_max 0.4101688 1.660005e-51
## CC.high_mean 0.3387289 1.149873e-34
## CC.mid_max 0.3125120 1.700196e-29
## CC.mid_mean 0.2523851 1.798179e-19
## WDir.900mb_mean 0.2387089 1.584530e-17
## Sunshine_sum 0.1905041 1.345970e-11
## Tot.Prcp_sum 0.1496478 1.195601e-07
## CC.low_min -0.1370961 1.261785e-06
## CC.low_mean -0.1520562 7.433425e-08
## Humidity_max -0.3903507 2.092512e-46
## Humidity_min -0.4485643 2.032901e-62
## Humidity_mean -0.5209515 3.327271e-87
## MSL_min -0.5487371 2.011197e-98
## MSL_mean -0.5921376 3.283794e-118
## MSL_max -0.6129165 8.518604e-129
##
## $Dim.3$quali
## R2 p.value
## pluie.demain 0.07172344 8.232673e-22
##
## $Dim.3$category
## Estimate p.value
## pluie.demain=TRUE 0.5606199 8.232673e-22
## pluie.demain=FALSE -0.5606199 8.232673e-22Graphique de corrélation avec corrplot
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## Temp_mean -0.40 0.07 0.78 0.26 0.23
## Humidity_mean 0.41 0.49 -0.52 0.06 0.01
## MSL_mean -0.41 -0.22 -0.59 0.34 0.19
## Tot.Prcp_sum 0.49 0.34 0.15 0.00 0.26
## Snow_sum 0.16 0.07 -0.10 -0.05 0.04
## Tot.CC_mean 0.74 0.54 -0.01 0.04 -0.02
## CC.high_mean 0.38 0.43 0.34 -0.25 0.05
## CC.mid_mean 0.66 0.48 0.25 -0.25 0.04
## CC.low_mean 0.69 0.47 -0.15 0.15 0.09
## Sunshine_sum -0.73 -0.47 0.19 0.08 0.06
## Radiation_sum -0.62 -0.21 0.53 0.18 0.04
## WS.10m_mean 0.84 -0.45 0.06 0.04 0.05
## WDir.10m_mean 0.39 0.26 0.05 0.75 -0.02
## WS.80m_mean 0.84 -0.46 0.05 0.00 0.05
## WDir.80m_mean 0.39 0.28 0.06 0.75 -0.03
## WS.900mb_mean 0.80 -0.41 0.02 -0.01 0.03
## WDir.900mb_mean 0.37 0.14 0.24 0.50 -0.04
## WGust_mean 0.86 -0.42 0.11 0.04 0.05
## Temp_max -0.47 0.01 0.76 0.18 0.20
## Temp_min -0.29 0.15 0.75 0.34 0.26
## Humidity_max 0.18 0.42 -0.39 -0.01 -0.10
## Humidity_min 0.44 0.45 -0.45 0.13 0.09
## MSL_max -0.33 -0.24 -0.61 0.31 0.18
## MSL_min -0.46 -0.18 -0.55 0.36 0.21
## Tot.CC_max 0.46 0.41 0.06 0.19 -0.42
## Tot.CC_min 0.39 0.45 -0.10 -0.12 0.64
## CC.high_max 0.32 0.42 0.41 -0.10 -0.14
## CC.high_min 0.18 0.14 0.08 -0.19 0.26
## CC.mid_max 0.50 0.40 0.31 -0.05 -0.27
## CC.mid_min 0.27 0.25 0.04 -0.23 0.46
## CC.low_max 0.52 0.38 -0.04 0.25 -0.34
## CC.low_min 0.27 0.36 -0.14 -0.07 0.62
## WS.10m_max 0.82 -0.37 0.09 0.01 -0.07
## WS.10m_min 0.68 -0.47 0.01 0.06 0.16
## WS.80m_max 0.80 -0.38 0.07 -0.08 -0.07
## WS.80m_min 0.69 -0.47 0.00 0.09 0.16
## WS.900mb_max 0.77 -0.32 0.05 -0.08 -0.09
## WS.900mb_min 0.69 -0.46 0.02 0.10 0.12
## WGust_max 0.81 -0.31 0.11 -0.03 -0.09
## WGust_min 0.71 -0.46 0.07 0.12 0.15
- Axe 1: On remarque que les 3 groupes de variables de type Vent: Wind Speed/Win dir et Win Gust. Et ont toutes une “corrélation” forte et >0 avec la 1ère dimension. De même pour Cloud Cover mais de moindre intensité. Opposition entre les variables de types Vents + Cloud Cover et les autres. Les autres types de variables: MSL/Sunshine duration/Radiation/Temperature ont toutes une “corrélation” <0.
- Axe 2: Les données de type Cloud Cover ont une “corrélation” >0 assez forte. Les données de type Win Speed et Wind Dir ont une “corrélation” <0 du même ordre.
- Axe 3: Les variables de type température ont la plus forte “corrélation” >0. Les variables de type MSL ont une “corrélation” <0 importante.
- Axe 4: Peu de variables ont une corrélation importantes. Excepté pour les variables de type Wind Dir qui ont une corrélation relativement importante.
- Axe 5: Les variables de type Cloud Cover: CC.Low.min Tot.CC.min
3.3.2 Indice de qualité de la représentation cos2
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## Temp_mean 0.16 0.00 0.61 0.07 0.05
## Humidity_mean 0.17 0.24 0.27 0.00 0.00
## MSL_mean 0.17 0.05 0.35 0.12 0.04
## Tot.Prcp_sum 0.24 0.12 0.02 0.00 0.07
## Snow_sum 0.03 0.01 0.01 0.00 0.00
## Tot.CC_mean 0.55 0.29 0.00 0.00 0.00
## CC.high_mean 0.15 0.19 0.11 0.06 0.00
## CC.mid_mean 0.43 0.23 0.06 0.06 0.00
## CC.low_mean 0.48 0.22 0.02 0.02 0.01
## Sunshine_sum 0.53 0.22 0.04 0.01 0.00
## Radiation_sum 0.39 0.04 0.28 0.03 0.00
## WS.10m_mean 0.71 0.20 0.00 0.00 0.00
## WDir.10m_mean 0.15 0.07 0.00 0.56 0.00
## WS.80m_mean 0.70 0.22 0.00 0.00 0.00
## WDir.80m_mean 0.15 0.08 0.00 0.56 0.00
## WS.900mb_mean 0.63 0.17 0.00 0.00 0.00
## WDir.900mb_mean 0.13 0.02 0.06 0.25 0.00
## WGust_mean 0.74 0.18 0.01 0.00 0.00
## Temp_max 0.23 0.00 0.58 0.03 0.04
## Temp_min 0.08 0.02 0.56 0.11 0.07
## Humidity_max 0.03 0.18 0.15 0.00 0.01
## Humidity_min 0.19 0.20 0.20 0.02 0.01
## MSL_max 0.11 0.06 0.38 0.10 0.03
## MSL_min 0.21 0.03 0.30 0.13 0.04
## Tot.CC_max 0.21 0.17 0.00 0.04 0.17
## Tot.CC_min 0.15 0.21 0.01 0.02 0.42
## CC.high_max 0.10 0.18 0.17 0.01 0.02
## CC.high_min 0.03 0.02 0.01 0.04 0.07
## CC.mid_max 0.25 0.16 0.10 0.00 0.08
## CC.mid_min 0.07 0.06 0.00 0.05 0.21
## CC.low_max 0.27 0.14 0.00 0.06 0.11
## CC.low_min 0.08 0.13 0.02 0.00 0.39
## WS.10m_max 0.67 0.14 0.01 0.00 0.00
## WS.10m_min 0.46 0.22 0.00 0.00 0.02
## WS.80m_max 0.64 0.14 0.01 0.01 0.00
## WS.80m_min 0.47 0.22 0.00 0.01 0.03
## WS.900mb_max 0.60 0.10 0.00 0.01 0.01
## WS.900mb_min 0.48 0.21 0.00 0.01 0.02
## WGust_max 0.66 0.10 0.01 0.00 0.01
## WGust_min 0.51 0.21 0.00 0.01 0.02
3.3.3 contribution des variables à la construction des axes
## Dim.1 Dim.2 Dim.3 Dim.4 Dim.5
## Temp_mean 1.22 0.09 13.82 2.71 2.73
## Humidity_mean 1.30 4.48 6.19 0.16 0.00
## MSL_mean 1.29 0.86 8.00 4.92 1.84
## Tot.Prcp_sum 1.84 2.13 0.51 0.00 3.51
## Snow_sum 0.20 0.10 0.23 0.10 0.08
## Tot.CC_mean 4.22 5.40 0.00 0.05 0.01
## CC.high_mean 1.14 3.41 2.62 2.66 0.10
## CC.mid_mean 3.30 4.16 1.45 2.60 0.10
## CC.low_mean 3.69 4.06 0.53 0.90 0.38
## Sunshine_sum 4.07 4.02 0.83 0.28 0.20
## Radiation_sum 2.99 0.79 6.50 1.37 0.10
## WS.10m_mean 5.48 3.76 0.08 0.08 0.11
## WDir.10m_mean 1.14 1.27 0.05 23.09 0.01
## WS.80m_mean 5.37 3.96 0.05 0.00 0.14
## WDir.80m_mean 1.15 1.45 0.10 23.49 0.04
## WS.900mb_mean 4.88 3.07 0.01 0.01 0.06
## WDir.900mb_mean 1.03 0.34 1.30 10.35 0.08
## WGust_mean 5.66 3.22 0.25 0.05 0.13
## Temp_max 1.73 0.00 13.30 1.28 2.11
## Temp_min 0.63 0.43 12.88 4.69 3.43
## Humidity_max 0.25 3.24 3.48 0.01 0.50
## Humidity_min 1.50 3.69 4.59 0.70 0.45
## MSL_max 0.86 1.03 8.57 3.97 1.66
## MSL_min 1.63 0.62 6.87 5.31 2.27
## Tot.CC_max 1.61 3.08 0.08 1.48 8.80
## Tot.CC_min 1.15 3.79 0.21 0.65 21.13
## CC.high_max 0.80 3.28 3.84 0.45 0.93
## CC.high_min 0.25 0.35 0.16 1.52 3.35
## CC.mid_max 1.95 2.87 2.23 0.10 3.82
## CC.mid_min 0.55 1.17 0.04 2.29 10.58
## CC.low_max 2.09 2.66 0.05 2.54 5.81
## CC.low_min 0.58 2.34 0.43 0.20 19.71
## WS.10m_max 5.13 2.52 0.18 0.00 0.23
## WS.10m_min 3.51 4.10 0.00 0.15 1.26
## WS.80m_max 4.93 2.64 0.11 0.24 0.24
## WS.80m_min 3.61 4.11 0.00 0.35 1.34
## WS.900mb_max 4.59 1.90 0.05 0.25 0.44
## WS.900mb_min 3.68 3.93 0.01 0.39 0.79
## WGust_max 5.05 1.78 0.28 0.04 0.42
## WGust_min 3.92 3.89 0.11 0.55 1.12
la contribution est très répartie entre les variables et diffuse pour les 2 premiers axes. On trouve surtout les variables de type Wind Speed et Wind Gust pour le 1ère axe. Les variables de type Cloud Cover pour le second mais aussi pour le 1ère axe. Les variables de type température pour le 3ème axe contribution moyenne. Pour le 4ème axe contribution forte de 2 variables Win Dir 10/80 m. Les variables Cloud Cover (Nébulosité) contribuent au 5ème axe.
3.3.4 Détail - Qualité et Contribution des variables aux axes
Détail - Qualité et Contribution des variables au 1er axe
which(res.pca.meteoTrain$var$contrib[,1]>4.5)## WS.10m_mean WS.80m_mean WS.900mb_mean WGust_mean WS.10m_max
## 12 14 16 18 33
## WS.80m_max WS.900mb_max WGust_max
## 35 37 39- 1ere Groupe de variables du type Vent avec 35% de la contribution WGust_mean, WS.10m_mean, WS.80m_mean, WS.10m_max, WGust_max, WS.80m_max, WS.900mb_mean, WS.900mb_max (>4.5%) plus de 50% pour les données de type Wind Speed / Wind Gust
- 2ème groupe de variables de type Cloud Coverrage avec plus de 10% pour les données de type Cloud Cover - Tot.CC_mean 5%,CC.low_mean 3%, CC.mid_mean 3%
- Autres: 7% Sunshine_sum 4%,Radiation_sum 3% mais en corrélation négative
Bonne représentation des variables de type Wind Speed et Wind Gust sur l’axe 1:
which(res.pca.meteoTrain$var$cos2[,1]>0.6)## WS.10m_mean WS.80m_mean WS.900mb_mean WGust_mean WS.10m_max
## 12 14 16 18 33
## WS.80m_max WGust_max
## 35 39On retiendra en premier lieu les variables: WS.10m_mean, WS.80m_mean, WS.10m_max, WGust_max, WS.80m_max, WS.900mb_mean, WS.900mb_max Puis: Tot.CC_mean, CC.low_mean, CC.mid_mean Et enfin: Sunshine_sum, Radiation_sum
Détail - Qualité et Contribution des variables au 2ème axe
Aucunes des variables ne se détachent nettement.
which(res.pca.meteoTrain$var$contrib[,2]>3)## Humidity_mean Tot.CC_mean CC.high_mean CC.mid_mean CC.low_mean
## 2 6 7 8 9
## Sunshine_sum WS.10m_mean WS.80m_mean WS.900mb_mean WGust_mean
## 10 12 14 16 18
## Humidity_max Humidity_min Tot.CC_max Tot.CC_min CC.high_max
## 21 22 25 26 27
## WS.10m_min WS.80m_min WS.900mb_min WGust_min
## 34 36 38 40Ces variables sont plutôt mal représentées sur l’axe 2 :
which(res.pca.meteoTrain$var$cos2[,2]>0.25 )## Tot.CC_mean
## 6which(res.pca.meteoTrain$var$cos2[,2]>0.2 )## Humidity_mean Tot.CC_mean CC.mid_mean CC.low_mean Sunshine_sum
## 2 6 8 9 10
## WS.10m_mean WS.80m_mean Humidity_min Tot.CC_min WS.10m_min
## 12 14 22 26 34
## WS.80m_min WS.900mb_min WGust_min
## 36 38 40On retiendra principalement: Tot.CC_mean + Humidity_mean + Humidity_min On pourrait ajouter: CC.mid_mean + CC.low_mean et WS.10m_min + WS.80m_min A noter qu’un nombre important de variables du 1er axe contribuent aussi au second. 
Détail - Qualité et Contribution des variables au 3ème axe
Temp_mean,Temp_max,Temp_min : contribuent à 40% à 50% à la construction de cet axe.
which(res.pca.meteoTrain$var$contrib[,3]>12)## Temp_mean Temp_max Temp_min
## 1 19 20Ces variables sont aussi les mieux représentés sur cette axe. A retenir pour l’explication de l’axe 3. Vient ensuite les variables MSL_max, MSL_mean, MSL_min qui expliquent +20% de la contribution. Ces variables de type MSL (Pression) sont assez bien représentées avec un cos2 > 0.3 On pourrait envisager de les retenir. Autres: Radiation_sum / Humidity_mean, Humidity_min / CC.high_max, CC.high_mean, CC.mid_max
Détail - Qualité et Contribution des variables au le 4ème axe
Cette dimension est expliquée par la direction du vent essentiellement WDir.80m_mean, WDir.10m_mean, WDir.900mb_mean 50% de la contribution
which(res.pca.meteoTrain$var$contrib[,4]>22)## WDir.10m_mean WDir.80m_mean
## 13 15Ensuite en 2ème groupe explicatif viennent les variables de type MSL : MSL_min, MSL_mean, MSL_max avec 15% de la contribution. Mais sont mals représentées avec un cos2 proche de 0.1 Autres: Temp_min / CC.low_max, CC.mid_min, Tot.CC_max Pour l’axe 4 comme variable explicative on va garder les variables de type Wind direction: WDir.80m_mean, WDir.10m_mean, WDir.900mb_mean
Détail - Qualité et Contribution des variables au le 5ème axe
Les variables Cloud Cover (Nébulosité) contribuent le plus au 5ème axe, à hauteur de 50%. Tot.CC_min + CC.mid_min + CC.low_min
which(res.pca.meteoTrain$var$contrib[,5]>10)## Tot.CC_min CC.mid_min CC.low_min
## 26 30 32
Remarque: (Sans l’option quanti.sup = [1:4] on aurait Year + X qui représente 60% de l’axe) C’est la différence essentielle.
3.4 Conclusion de l’analyse PCA - Choix des variables et modèle retenu
3.4.1 Variables retenues à l’issue de l’analyse PCA
Pour expliquer les axes factoriels on va retenir:
- 1er axe : WGust_mean + WS.10m_mean + WS.80m_mean + WS.10m_max + WGust_max + WS.80m_max + WS.900mb_mean + WS.900mb_max Tot.CC_mean, CC.low_mean, CC.mid_mean Sunshine_sum, Radiation_sum
- 2ème axe : Tot.CC_mean + Humidity_mean On pourrait ajouter: CC.mid_mean + CC.low_mean / Humidity_min / WS.10m_min + WS.80m_min A noter qu’un nombre important de variables du 1er axe contribuent aussi au second.
- 3ème axe : Temp_mean + Temp_max + Temp_min
- 4ème axe : WDir.80m_mean + WDir.10m_mean + WDir.900mb_mean
- 5ème axe : Tot.CC_min + CC.mid_min + CC.low_min (Sans l’option quali.sup = [1:4] on aurait Year + X qui représente 60% de l’axe)
On pourrait ajouter aussi les variables suivantes, qui peuvent apporter une certaine diversification. Et qui apparaissent en 2ème ordre dans l’analyse: MSL_min + MSL_max + MSL_mean et Sunshine_sum Ainsi que CC.mid_max, Tot.Prcp_sum, Snow_sum et Radiation_sum pour diversification.
3.4.2 Modèle issue de la sélection de variables à partir de l’analyse PCA
# Avec l'option quali.sup=[1:4]
m_PCA = glm(pluie.demain ~
WGust_mean + WGust_max #1er axe
+ WS.10m_mean + WS.80m_mean + WS.900mb_mean #1er axe
+ WS.10m_max + WS.80m_max + WS.900mb_max #1er axe
+ Tot.CC_mean + Humidity_mean + Humidity_min + #2eme axe
+ WS.10m_min + WS.80m_min + #2eme axe
+ Temp_mean + Temp_max + Temp_min + #3eme axe
+ WDir.80m_mean + WDir.10m_mean + WDir.900mb_mean + #4eme axe
+ Tot.CC_min + CC.mid_min + CC.low_min #5eme axe
+ CC.mid_max + MSL_max + MSL_mean + MSL_min #divers
+ Sunshine_sum + Tot.Prcp_sum + Snow_sum + Radiation_sum #divers
,family = binomial, data = meteoTrain )
summary(m_PCA)##
## Call:
## glm(formula = pluie.demain ~ WGust_mean + WGust_max + WS.10m_mean +
## WS.80m_mean + WS.900mb_mean + WS.10m_max + WS.80m_max + WS.900mb_max +
## Tot.CC_mean + Humidity_mean + Humidity_min + +WS.10m_min +
## WS.80m_min + +Temp_mean + Temp_max + Temp_min + +WDir.80m_mean +
## WDir.10m_mean + WDir.900mb_mean + +Tot.CC_min + CC.mid_min +
## CC.low_min + CC.mid_max + MSL_max + MSL_mean + MSL_min +
## Sunshine_sum + Tot.Prcp_sum + Snow_sum + Radiation_sum, family = binomial,
## data = meteoTrain)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.4848 -0.8890 -0.2178 0.9211 2.3324
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 8.483e+01 1.174e+01 7.227 4.94e-13 ***
## WGust_mean 3.122e-02 2.951e-02 1.058 0.29011
## WGust_max 1.224e-02 1.485e-02 0.825 0.40961
## WS.10m_mean -5.978e-02 9.352e-02 -0.639 0.52267
## WS.80m_mean -9.866e-02 6.442e-02 -1.532 0.12560
## WS.900mb_mean -3.314e-03 1.519e-02 -0.218 0.82735
## WS.10m_max 2.713e-02 3.238e-02 0.838 0.40217
## WS.80m_max 1.185e-03 2.674e-02 0.044 0.96464
## WS.900mb_max 7.873e-03 9.487e-03 0.830 0.40663
## Tot.CC_mean 2.749e-03 7.318e-03 0.376 0.70719
## Humidity_mean -1.029e-02 1.792e-02 -0.574 0.56564
## Humidity_min 1.473e-02 1.411e-02 1.044 0.29642
## WS.10m_min 1.074e-01 5.448e-02 1.971 0.04872 *
## WS.80m_min -6.750e-04 3.874e-02 -0.017 0.98610
## Temp_mean 1.883e-01 1.546e-01 1.218 0.22311
## Temp_max 2.224e-02 9.098e-02 0.244 0.80686
## Temp_min -2.232e-01 8.101e-02 -2.755 0.00587 **
## WDir.80m_mean -8.179e-03 5.200e-03 -1.573 0.11570
## WDir.10m_mean 6.514e-03 5.016e-03 1.299 0.19402
## WDir.900mb_mean 6.206e-03 1.350e-03 4.596 4.30e-06 ***
## Tot.CC_min 2.009e-02 8.295e-03 2.422 0.01543 *
## CC.mid_min 6.779e-03 1.059e-02 0.640 0.52218
## CC.low_min -2.200e-02 8.402e-03 -2.619 0.00882 **
## CC.mid_max 9.541e-03 2.158e-03 4.421 9.81e-06 ***
## MSL_max -1.647e-01 7.249e-02 -2.272 0.02307 *
## MSL_mean 2.033e-01 1.280e-01 1.588 0.11227
## MSL_min -1.246e-01 6.569e-02 -1.898 0.05776 .
## Sunshine_sum -4.711e-04 8.269e-04 -0.570 0.56884
## Tot.Prcp_sum 2.971e-02 2.246e-02 1.323 0.18580
## Snow_sum 8.477e-01 4.795e-01 1.768 0.07710 .
## Radiation_sum 2.922e-05 8.693e-05 0.336 0.73680
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1719.0 on 1239 degrees of freedom
## Residual deviance: 1354.7 on 1209 degrees of freedom
## AIC: 1416.7
##
## Number of Fisher Scoring iterations: 6anova(m_PCA,test="LRT")## Analysis of Deviance Table
##
## Model: binomial, link: logit
##
## Response: pluie.demain
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL 1239 1719.0
## WGust_mean 1 47.848 1238 1671.1 4.605e-12 ***
## WGust_max 1 24.881 1237 1646.2 6.097e-07 ***
## WS.10m_mean 1 1.172 1236 1645.1 0.2790756
## WS.80m_mean 1 6.824 1235 1638.2 0.0089958 **
## WS.900mb_mean 1 1.731 1234 1636.5 0.1882495
## WS.10m_max 1 0.146 1233 1636.4 0.7025441
## WS.80m_max 1 1.117 1232 1635.3 0.2905674
## WS.900mb_max 1 9.332 1231 1625.9 0.0022523 **
## Tot.CC_mean 1 58.243 1230 1567.7 2.317e-14 ***
## Humidity_mean 1 11.336 1229 1556.3 0.0007604 ***
## Humidity_min 1 0.379 1228 1556.0 0.5379390
## WS.10m_min 1 6.105 1227 1549.9 0.0134775 *
## WS.80m_min 1 0.534 1226 1549.3 0.4648743
## Temp_mean 1 18.488 1225 1530.8 1.710e-05 ***
## Temp_max 1 6.125 1224 1524.7 0.0133269 *
## Temp_min 1 9.443 1223 1515.3 0.0021195 **
## WDir.80m_mean 1 1.226 1222 1514.0 0.2681944
## WDir.10m_mean 1 1.396 1221 1512.7 0.2373889
## WDir.900mb_mean 1 34.981 1220 1477.7 3.329e-09 ***
## Tot.CC_min 1 1.493 1219 1476.2 0.2218164
## CC.mid_min 1 3.463 1218 1472.7 0.0627396 .
## CC.low_min 1 6.523 1217 1466.2 0.0106477 *
## CC.mid_max 1 35.695 1216 1430.5 2.308e-09 ***
## MSL_max 1 65.385 1215 1365.1 6.161e-16 ***
## MSL_mean 1 0.273 1214 1364.8 0.6016202
## MSL_min 1 3.182 1213 1361.7 0.0744682 .
## Sunshine_sum 1 0.315 1212 1361.3 0.5748339
## Tot.Prcp_sum 1 1.823 1211 1359.5 0.1769138
## Snow_sum 1 4.666 1210 1354.8 0.0307574 *
## Radiation_sum 1 0.113 1209 1354.7 0.7367225
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 13.4.3 Corrélations entre variables explicatives
La première dimension qui représente plus de 32% est expliqué essentiellement par des variables de catégorie vents. On peut essayer de regarder les corrélations 2/2 des groupes de variables sélectionnées avec PCA.
Peut être est-il intéressant de rechercher des corrélations au sein de certains groupes de variables et d’en éliminer certaines. Les variables de type vents sont plutôt très corrélées. En particulier WS.80m et WS10m Max et Mean Par la suite, lors de la selection de modèles on verra que ces variables de type vent-vitesse sont éliminées en grande partie. On n’en retient généralement 2 ou 3, voire une seule.
4. Sélection de modèles - Choix de modèles (méthode exhaustive)
Face à la grande quantité de variables explicatives, on cherche à les sélectionner de manière automatique.
L’utilisation de la méthode regsubset du package leaps à partir du modèle saturé s’est montrée infructueuse. La convergence exacte des critères Cp /BIC/IC… est succpecte. Cette méthode et les modèles obtenus ne sont pas satifaisants. Y-a t’il trop de variables pour que la méthode regsubset soit significative? Ou bien est-ce cette méthode n’est pas adaptée aux modèles linéaire généralisés tel que logit. En effet on n’a étudié cette méthode, que dans le cadre des modèles linéaire gaussien.
Comme alternative on va utiliser la méthode pas a pas: step du package MASS. Qui est basé sur une minimsation du critère AIC. On commence par la méthode progressive, puis descenadante et ascenadante à parir du modèle saturé.
Dans un deuxième temps on utilisera à nouveaux ces méthodes mais à partir du modèle réduit issue de l’analyse PCA: m_PCA. En effet beaucoup de méthodes de sélection sont applicable à un nombre de variable réduit de l’ordre de 30 ou moins. On ne peut donc les utiliser à partir du modèle saturé à 44 variables. On essaiera aussi à partir de la fonction leaps d’utiliser les critères de R2aj et le Cp de Mallow.
4.1 Méthodes step du package MASS à partir du modèle saturé
On reprend le modèle saturé obtenu au §2: m_sature
m_sature = glm(formula = pluie.demain ~ (X+Year+Month+Day+
Temp_mean + Humidity_mean + MSL_mean + Tot.Prcp_sum + Snow_sum +
Tot.CC_mean + CC.high_mean + CC.mid_mean + CC.low_mean +
Sunshine_sum + Radiation_sum + WS.10m_mean + WDir.10m_mean +
WS.80m_mean + WDir.80m_mean + WS.900mb_mean + WDir.900mb_mean +
WGust_mean + Temp_max + Temp_min + Humidity_max + Humidity_min +
MSL_max + MSL_min + Tot.CC_max + Tot.CC_min + CC.high_max +
CC.high_min + CC.mid_max + CC.mid_min + CC.low_max + CC.low_min +
WS.10m_max + WS.10m_min + WS.80m_max + WS.80m_min + WS.900mb_max +
WS.900mb_min + WGust_max + WGust_min), family = binomial(link="logit"), data = meteoTrain)4.1.1 Méthode progressive - step backward-forward à partir de m_sature
library(MASS)
modele_step_BwdFwd <- step(m_sature, data=meteoTrain, direction="both")Modèle obtenu:
m_StepBwdFwd<-glm(formula = pluie.demain ~ Temp_mean + Snow_sum + CC.low_mean + WS.80m_mean +
WDir.900mb_mean + Temp_min + MSL_max + Tot.CC_max + Tot.CC_min +
CC.high_max + CC.mid_max + CC.low_min + WS.10m_min + WGust_max +
WS.900mb_max, family = binomial, data = meteoTrain)
AIC(m_StepBwdFwd)## [1] 1392.3724.1.2 Méthode descendante - step backward à partir de m_sature
modele_step_Bwd <- step(m_sature, data=meteoTrain, direction="backward")Modèle obtenu:
m_StepBwd<-glm(pluie.demain ~ Temp_mean + Snow_sum + CC.low_mean + WS.80m_mean +
WDir.900mb_mean + Temp_min + MSL_max + MSL_min + Tot.CC_min +
CC.high_max + CC.mid_max + CC.low_min + WS.10m_min + WGust_max,
data = meteoTrain, family = binomial)
#summary(m_StepBwd)
AIC(m_StepBwd)## [1] 1393.1644.1.3 Méthode ascendante - step forward à partir de m_sature
step(glm(pluie.demain ~1,data=meteoTrain, family = binomial), pluie.demain ~
X+Year+Month+Day+ Temp_mean + Humidity_mean + MSL_mean + Tot.Prcp_sum + Snow_sum +
Tot.CC_mean + CC.high_mean + CC.mid_mean + CC.low_mean +
Sunshine_sum + Radiation_sum + WS.10m_mean + WDir.10m_mean +
WS.80m_mean + WDir.80m_mean + WS.900mb_mean + WDir.900mb_mean +
WGust_mean + Temp_max + Temp_min + Humidity_max + Humidity_min +
MSL_max + MSL_min + Tot.CC_max + Tot.CC_min + CC.high_max +
CC.high_min + CC.mid_max + CC.mid_min + CC.low_max + CC.low_min +
WS.10m_max + WS.10m_min + WS.80m_max + WS.80m_min + WS.900mb_max +
WS.900mb_min + WGust_max + WGust_min,
data=meteoTrain, direction="forward")Modèle obtenu:
m_StepFwd<-glm(formula = pluie.demain ~ MSL_min + CC.mid_max + WDir.900mb_mean +
CC.high_mean + Snow_sum + MSL_max + MSL_mean + WS.900mb_max +
WS.80m_mean + WS.10m_min + WGust_max + Tot.CC_max, family = binomial,
data = meteoTrain)
#summary(m_StepFwd)
AIC(m_StepFwd)## [1] 1413.5884.2 Méthodes step du package MASS à partir du modèle pré-sélectionné avec l’analyse PCA
Reprenons le modèle PCA: m_PCA
m_PCA<-glm(pluie.demain ~ 1 +
WGust_mean + WGust_max #1er axe
+ WS.10m_mean + WS.80m_mean + WS.900mb_mean #1er axe
+ WS.10m_max + WS.80m_max + WS.900mb_max #1er axe
+ Tot.CC_mean + Humidity_mean + Humidity_min + #2eme axe
+ WS.10m_min + WS.80m_min + #2eme axe
+ Temp_mean + Temp_max + Temp_min + #3eme axe
+ WDir.80m_mean + WDir.10m_mean + WDir.900mb_mean + #4eme axe
+ Tot.CC_min + CC.low_min + CC.mid_min #5eme axe
+ CC.mid_max + MSL_max + MSL_mean + MSL_min #divers
+ Sunshine_sum + Tot.Prcp_sum + Snow_sum + Radiation_sum, #divers
family = binomial, data = meteoTrain )4.2.1 Méthode progressive - step backward-forward à partir de m_PCA
library(MASS)
modele_step_BwdFwdPCA <- step(m_PCA, data=meteoTrain, direction="both")Modèle obtenu:
m_StepBwdFwd_PCA<-glm(formula = pluie.demain ~ 1 + WGust_max + WS.80m_mean + WS.10m_min + Temp_mean +
Temp_min + WDir.900mb_mean + Tot.CC_min + CC.low_min + CC.mid_max +
MSL_max + MSL_mean + MSL_min + Sunshine_sum + Tot.Prcp_sum + Snow_sum,
family = binomial, data = meteoTrain)
AIC(m_StepBwdFwd_PCA)## [1] 1394.18BIC(m_StepBwdFwd_PCA)## [1] 1476.1464.2.2 Méthode descendante - step backward à partir de m_PCA
modele_step_BwdPCA <- step(m_PCA, data=meteoTrain, direction="backward")Modèle obtenu:
m_StepBwd_PCA<-glm(pluie.demain ~ WGust_max + WS.80m_mean + WS.10m_min + Temp_mean +
Temp_min + WDir.900mb_mean + Tot.CC_min + CC.low_min + CC.mid_max +
MSL_max + MSL_mean + MSL_min + Sunshine_sum + Tot.Prcp_sum + Snow_sum, data = meteoTrain, family = binomial)
AIC(m_StepBwd_PCA)## [1] 1394.18On obtient le même modèle qu’avec la procédure progressive: m_StepBwd_PCA = m_StepBwdFwd_PCA
4.2.3 Méthode ascendante - step forward à partir de m_PCA
step(glm(pluie.demain ~1,data=meteoTrain, family = binomial),
pluie.demain ~ 1 +
WGust_mean + WGust_max #1er axe
+ WS.10m_mean + WS.80m_mean + WS.900mb_mean #1er axe
+ WS.10m_max + WS.80m_max + WS.900mb_max #1er axe
+ Tot.CC_mean + Humidity_mean + Humidity_min + #2eme axe
+ WS.10m_min + WS.80m_min + #2eme axe
+ Temp_mean + Temp_max + Temp_min + #3eme axe
+ WDir.80m_mean + WDir.10m_mean + WDir.900mb_mean + #4eme axe
+ Tot.CC_min + CC.low_min + CC.mid_min #5eme axe
+ CC.mid_max + MSL_max + MSL_mean + MSL_min #divers
+ Sunshine_sum + Tot.Prcp_sum + Snow_sum + Radiation_sum, #divers
family = binomial, data=meteoTrain, direction="forward")Modèle obtenu:
m_StepFwd_PCA<-glm(formula = pluie.demain ~ MSL_min + CC.mid_max + WDir.900mb_mean +
Snow_sum + CC.mid_min + MSL_max + MSL_mean + WS.900mb_max +
WS.80m_mean + WS.10m_min + WGust_max + WDir.80m_mean,
family = binomial, data = meteoTrain)
AIC(m_StepFwd_PCA)## [1] 1416.7524.3 Autres méthodes: utilisation du package leaps
Pour réaliser cette sélection de variables on va utiliser la méthode leaps. Cette méthode n’est utilsable que pour au plus 31 variables. On est donc obligé de faire une première pré-selection. Pour celà on va partir des variables obtenues à partir de l’analyse PCA. Une alternative serait d’effectuer une pré-selection avec la méthode step (AIC) et poursuivre avec les autres critères.
Préparation:
Chargement de la librairie leaps et tableau des numéro de colonnes
4.3.1 Utilisation du R2Ajusté à partir de la méthode leaps
Covariables retenues:
print(names(meteoTrain)[colnum][best.mod])## [1] "Temp_mean" "MSL_mean" "Tot.Prcp_sum"
## [4] "Snow_sum" "Sunshine_sum" "WDir.10m_mean"
## [7] "WS.80m_mean" "WDir.80m_mean" "WDir.900mb_mean"
## [10] "Temp_min" "MSL_max" "MSL_min"
## [13] "Tot.CC_min" "CC.mid_max" "CC.low_min"
## [16] "WS.10m_min" "WS.900mb_max" "WGust_max"Modèle obtenu:
m_R2aj_PCA<-glm(pluie.demain ~
Temp_mean + MSL_mean + Tot.Prcp_sum + Snow_sum + Sunshine_sum +
WDir.10m_mean + WS.80m_mean + WDir.80m_mean + WDir.900mb_mean +
Temp_min + MSL_max + MSL_min + Tot.CC_min + CC.mid_max + CC.low_min +
WS.10m_min + WS.900mb_max + WGust_max, family = binomial, data = meteoTrain )
AIC(m_R2aj_PCA)## [1] 1396.347Le résultat du critère AIC est correcte.
4.3.2 Utilisation du Cp de Mallow à partir de la méthode leaps
Covariables retenues:
## [1] "Temp_mean" "MSL_mean" "Tot.Prcp_sum"
## [4] "Sunshine_sum" "WS.80m_mean" "WDir.900mb_mean"
## [7] "Temp_min" "MSL_max" "MSL_min"
## [10] "Tot.CC_min" "CC.mid_max" "CC.low_min"
## [13] "WS.10m_min" "WGust_max"Modèle obtenu:
m_Cp_PCA<-glm(pluie.demain ~
Temp_mean + MSL_mean + Tot.Prcp_sum + Sunshine_sum + WS.80m_mean + WDir.900mb_mean
+ Temp_min + MSL_max + MSL_min + Tot.CC_min + CC.mid_max + CC.low_min
+ WGust_max + WS.10m_min ,family = binomial, data = meteoTrain )
AIC(m_Cp_PCA)## [1] 1396.6794.3.3 Tentative d’utilisation de regsubsets à partir du modèle réduit m_PCA
Là encore le résultat n’est pas satisfaisant. 
4.4 Autres procédures automatiques adaptées aux modèles glm
4.4.1 Le package bestglm
Apparamment limité à 15 variables explicatives. La solution serait de réduire le nombre de variables à 15 par la méthode step. Puis d’utiliser la méthode bestglm.
http://www2.uaem.mx/r-mirror/web/packages/bestglm/vignettes/bestglm.pdf
4.4.2 Les packages glmnet et BeSS
Sortent du cadre du cours.
4.4.3 Le package glmulti
Limité a 32 variables explicatives. On va pouvoir l’utiliser à partir des variables préselectionnées lors de l’analyse PCA. Le code est désactivé (en commentaire) car la méthode est très coûteuse en temps.
Le coût en temps étant vraiment trop important, il a fallut éliminer encore des variables. Les modèles obtenus ont un moins bon AIC que nos autres modèles.
Modèles obtenus à partir du package glmulti:
m_glmulti_PCA1<-glm(formula = pluie.demain~
WS.80m_mean + WS.900mb_max + WGust_max +
Temp_mean + Temp_min +
WDir.900mb_mean +
MSL_max +
Tot.CC_mean + CC.mid_max + Tot.CC_min + CC.low_min, family = binomial,
data = meteoTrain)
AIC(m_glmulti_PCA1)## [1] 1405.837BIC(m_glmulti_PCA1)## [1] 1467.312m_glmulti_PCA2<-glm(formula = pluie.demain~
WS.10m_mean + WS.80m_mean + WGust_max +
Temp_mean + Temp_min +
WDir.80m_mean + WDir.10m_mean + WDir.900mb_mean +
MSL_max + MSL_min + Sunshine_sum +
CC.mid_max + Tot.CC_min + CC.low_min, family = binomial,
data = meteoTrain)
AIC(m_glmulti_PCA2)## [1] 1407.505BIC(m_glmulti_PCA2)## [1] 1484.348m_glmulti_PCA3<-glm(formula = pluie.demain~
WS.10m_mean + WS.80m_mean + WS.900mb_max +
Temp_mean + Temp_min +
WDir.900mb_mean +
MSL_max +
Tot.CC_mean + CC.mid_max + Tot.CC_min + CC.low_min, family = binomial,
data = meteoTrain)
AIC(m_glmulti_PCA3)## [1] 1407.399BIC(m_glmulti_PCA3)## [1] 1468.873Vraissemblablement la préselection s’est révélée trop restrictive. On ne va donc pas garder ces modèles.
4.5 Résumé des caractéristiques des modèles retenus
4.5.1 Les modèles obtenus à partir du modèle saturé
Le modèle saturé ou complet
#Modèle saturé
m_sature<-glm(formula = pluie.demain ~ (X+Year+Month+Day+
Temp_mean + Humidity_mean + MSL_mean + Tot.Prcp_sum + Snow_sum +
Tot.CC_mean + CC.high_mean + CC.mid_mean + CC.low_mean +
Sunshine_sum + Radiation_sum + WS.10m_mean + WDir.10m_mean +
WS.80m_mean + WDir.80m_mean + WS.900mb_mean + WDir.900mb_mean +
WGust_mean + Temp_max + Temp_min + Humidity_max + Humidity_min +
MSL_max + MSL_min + Tot.CC_max + Tot.CC_min + CC.high_max +
CC.high_min + CC.mid_max + CC.mid_min + CC.low_max + CC.low_min +
WS.10m_max + WS.10m_min + WS.80m_max + WS.80m_min + WS.900mb_max +
WS.900mb_min + WGust_max + WGust_min), family = binomial, data = meteoTrain)
summary(m_sature)##
## Call:
## glm(formula = pluie.demain ~ (X + Year + Month + Day + Temp_mean +
## Humidity_mean + MSL_mean + Tot.Prcp_sum + Snow_sum + Tot.CC_mean +
## CC.high_mean + CC.mid_mean + CC.low_mean + Sunshine_sum +
## Radiation_sum + WS.10m_mean + WDir.10m_mean + WS.80m_mean +
## WDir.80m_mean + WS.900mb_mean + WDir.900mb_mean + WGust_mean +
## Temp_max + Temp_min + Humidity_max + Humidity_min + MSL_max +
## MSL_min + Tot.CC_max + Tot.CC_min + CC.high_max + CC.high_min +
## CC.mid_max + CC.mid_min + CC.low_max + CC.low_min + WS.10m_max +
## WS.10m_min + WS.80m_max + WS.80m_min + WS.900mb_max + WS.900mb_min +
## WGust_max + WGust_min), family = binomial, data = meteoTrain)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.5638 -0.8784 -0.1965 0.9022 2.5003
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -7.342e+03 8.945e+04 -0.082 0.934580
## X -1.000e-02 1.218e-01 -0.082 0.934543
## Year 3.696e+00 4.449e+01 0.083 0.933794
## Month 3.005e-01 3.705e+00 0.081 0.935348
## Day 2.857e-03 1.218e-01 0.023 0.981280
## Temp_mean 1.603e-01 1.625e-01 0.987 0.323670
## Humidity_mean -3.515e-02 3.002e-02 -1.171 0.241580
## MSL_mean 2.138e-01 1.319e-01 1.620 0.105151
## Tot.Prcp_sum 2.436e-02 2.553e-02 0.954 0.339952
## Snow_sum 8.294e-01 4.838e-01 1.714 0.086442 .
## Tot.CC_mean -2.058e-03 1.208e-02 -0.170 0.864760
## CC.high_mean 4.740e-03 6.811e-03 0.696 0.486416
## CC.mid_mean -4.237e-03 6.529e-03 -0.649 0.516402
## CC.low_mean 7.270e-03 8.276e-03 0.878 0.379684
## Sunshine_sum -6.158e-04 8.673e-04 -0.710 0.477738
## Radiation_sum 4.286e-05 9.846e-05 0.435 0.663331
## WS.10m_mean -6.170e-02 9.611e-02 -0.642 0.520907
## WDir.10m_mean 7.129e-03 5.086e-03 1.402 0.161038
## WS.80m_mean -1.122e-01 6.587e-02 -1.703 0.088580 .
## WDir.80m_mean -8.893e-03 5.267e-03 -1.688 0.091318 .
## WS.900mb_mean 1.770e-02 2.291e-02 0.772 0.439835
## WDir.900mb_mean 6.630e-03 1.378e-03 4.812 1.49e-06 ***
## WGust_mean 3.863e-02 3.494e-02 1.106 0.268836
## Temp_max 3.185e-02 9.420e-02 0.338 0.735301
## Temp_min -2.117e-01 8.540e-02 -2.479 0.013177 *
## Humidity_max 2.138e-02 1.945e-02 1.099 0.271635
## Humidity_min 2.470e-02 1.649e-02 1.498 0.134246
## MSL_max -1.794e-01 7.451e-02 -2.407 0.016075 *
## MSL_min -1.245e-01 6.741e-02 -1.847 0.064818 .
## Tot.CC_max -5.009e-03 4.728e-03 -1.059 0.289431
## Tot.CC_min 1.908e-02 8.636e-03 2.209 0.027179 *
## CC.high_max 2.626e-03 2.754e-03 0.954 0.340300
## CC.high_min -3.508e-02 2.763e-02 -1.270 0.204157
## CC.mid_max 1.106e-02 3.088e-03 3.583 0.000339 ***
## CC.mid_min 1.017e-02 1.149e-02 0.885 0.375956
## CC.low_max -4.016e-04 3.323e-03 -0.121 0.903795
## CC.low_min -2.214e-02 8.749e-03 -2.530 0.011399 *
## WS.10m_max 2.591e-02 3.318e-02 0.781 0.434926
## WS.10m_min 1.256e-01 5.829e-02 2.155 0.031139 *
## WS.80m_max 2.154e-03 2.723e-02 0.079 0.936942
## WS.80m_min 1.233e-02 4.055e-02 0.304 0.761117
## WS.900mb_max 1.423e-03 1.079e-02 0.132 0.895088
## WS.900mb_min -1.286e-02 1.720e-02 -0.748 0.454501
## WGust_max 1.196e-02 1.705e-02 0.702 0.482980
## WGust_min -2.189e-02 2.625e-02 -0.834 0.404328
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1719.0 on 1239 degrees of freedom
## Residual deviance: 1341.8 on 1195 degrees of freedom
## AIC: 1431.8
##
## Number of Fisher Scoring iterations: 6Les 3 modèles obtenus par minimsation du critère AIC de la méthode step - (Both, Backward et Forward)
(A patir du modèle saturé m_sature)
#Modèle step backward/forward: m_StepBwdFwd à partir de m_sature
m_StepBwdFwd<-glm(formula = pluie.demain ~ Temp_mean + Snow_sum + CC.low_mean + WS.80m_mean +
WDir.900mb_mean + Temp_min + MSL_max + Tot.CC_max + Tot.CC_min +
CC.high_max + CC.mid_max + CC.low_min + WS.10m_min + WGust_max +
WS.900mb_max, family = binomial, data = meteoTrain)
summary(m_StepBwdFwd)##
## Call:
## glm(formula = pluie.demain ~ Temp_mean + Snow_sum + CC.low_mean +
## WS.80m_mean + WDir.900mb_mean + Temp_min + MSL_max + Tot.CC_max +
## Tot.CC_min + CC.high_max + CC.mid_max + CC.low_min + WS.10m_min +
## WGust_max + WS.900mb_max, family = binomial, data = meteoTrain)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.5216 -0.8908 -0.2197 0.9103 2.4453
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 85.602525 11.344689 7.546 4.50e-14 ***
## Temp_mean 0.219313 0.053123 4.128 3.65e-05 ***
## Snow_sum 0.823033 0.453487 1.815 0.069540 .
## CC.low_mean 0.008459 0.003489 2.425 0.015321 *
## WS.80m_mean -0.097198 0.022211 -4.376 1.21e-05 ***
## WDir.900mb_mean 0.005134 0.001071 4.795 1.63e-06 ***
## Temp_min -0.236706 0.056121 -4.218 2.47e-05 ***
## MSL_max -0.086822 0.011023 -7.876 3.37e-15 ***
## Tot.CC_max -0.005083 0.003423 -1.485 0.137476
## Tot.CC_min 0.020511 0.007537 2.721 0.006504 **
## CC.high_max 0.004111 0.002102 1.956 0.050475 .
## CC.mid_max 0.010046 0.002636 3.810 0.000139 ***
## CC.low_min -0.022423 0.008075 -2.777 0.005490 **
## WS.10m_min 0.094931 0.029302 3.240 0.001196 **
## WGust_max 0.024525 0.010186 2.408 0.016055 *
## WS.900mb_max 0.010142 0.005932 1.710 0.087336 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1719.0 on 1239 degrees of freedom
## Residual deviance: 1360.4 on 1224 degrees of freedom
## AIC: 1392.4
##
## Number of Fisher Scoring iterations: 6#Modèle step backward: m_StepBwd à partir de m_sature
m_StepBwd<-glm(pluie.demain ~ Temp_mean + Snow_sum + CC.low_mean + WS.80m_mean +
WDir.900mb_mean + Temp_min + MSL_max + MSL_min + Tot.CC_min +
CC.high_max + CC.mid_max + CC.low_min + WS.10m_min + WGust_max,
data = meteoTrain, family = binomial)
summary(m_StepBwd)##
## Call:
## glm(formula = pluie.demain ~ Temp_mean + Snow_sum + CC.low_mean +
## WS.80m_mean + WDir.900mb_mean + Temp_min + MSL_max + MSL_min +
## Tot.CC_min + CC.high_max + CC.mid_max + CC.low_min + WS.10m_min +
## WGust_max, family = binomial, data = meteoTrain)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.4573 -0.9021 -0.2122 0.9101 2.3823
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 85.090520 11.287800 7.538 4.76e-14 ***
## Temp_mean 0.216246 0.053025 4.078 4.54e-05 ***
## Snow_sum 0.823612 0.453345 1.817 0.069256 .
## CC.low_mean 0.006911 0.003370 2.050 0.040329 *
## WS.80m_mean -0.091502 0.021894 -4.179 2.92e-05 ***
## WDir.900mb_mean 0.005082 0.001067 4.764 1.90e-06 ***
## Temp_min -0.230447 0.056343 -4.090 4.31e-05 ***
## MSL_max -0.050988 0.025439 -2.004 0.045030 *
## MSL_min -0.035662 0.023854 -1.495 0.134909
## Tot.CC_min 0.022027 0.007507 2.934 0.003343 **
## CC.high_max 0.004023 0.002086 1.928 0.053794 .
## CC.mid_max 0.008237 0.002358 3.493 0.000478 ***
## CC.low_min -0.023298 0.008058 -2.891 0.003838 **
## WS.10m_min 0.097444 0.029338 3.321 0.000896 ***
## WGust_max 0.031579 0.008583 3.679 0.000234 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1719.0 on 1239 degrees of freedom
## Residual deviance: 1363.2 on 1225 degrees of freedom
## AIC: 1393.2
##
## Number of Fisher Scoring iterations: 6#Modèle step forward: m_StepFwd à partir de m_sature
m_StepFwd<-glm(formula = pluie.demain ~ MSL_min + CC.mid_max + WDir.900mb_mean +
CC.high_mean + Snow_sum + MSL_max + MSL_mean + WS.900mb_max +
WS.80m_mean + WS.10m_min + WGust_max + Tot.CC_max, family = binomial,
data = meteoTrain)
summary(m_StepFwd)##
## Call:
## glm(formula = pluie.demain ~ MSL_min + CC.mid_max + WDir.900mb_mean +
## CC.high_mean + Snow_sum + MSL_max + MSL_mean + WS.900mb_max +
## WS.80m_mean + WS.10m_min + WGust_max + Tot.CC_max, family = binomial,
## data = meteoTrain)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.3781 -0.9143 -0.2713 0.9443 2.4202
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 91.6383029 10.7870477 8.495 < 2e-16 ***
## MSL_min -0.1048518 0.0616240 -1.701 0.08885 .
## CC.mid_max 0.0111706 0.0024313 4.594 4.34e-06 ***
## WDir.900mb_mean 0.0039716 0.0009649 4.116 3.85e-05 ***
## CC.high_mean 0.0093876 0.0037377 2.512 0.01202 *
## Snow_sum 0.8854241 0.4492552 1.971 0.04874 *
## MSL_max -0.1308116 0.0676997 -1.932 0.05333 .
## MSL_mean 0.1441538 0.1200402 1.201 0.22980
## WS.900mb_max 0.0070727 0.0058742 1.204 0.22858
## WS.80m_mean -0.0874552 0.0219211 -3.990 6.62e-05 ***
## WS.10m_min 0.0829775 0.0288147 2.880 0.00398 **
## WGust_max 0.0242992 0.0100767 2.411 0.01589 *
## Tot.CC_max -0.0050841 0.0031582 -1.610 0.10744
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1719.0 on 1239 degrees of freedom
## Residual deviance: 1387.6 on 1227 degrees of freedom
## AIC: 1413.6
##
## Number of Fisher Scoring iterations: 6Comparaison avec le modèle “saturé” d’origine:
anova(m_sature,m_StepBwdFwd,test="LRT")
anova(m_sature,m_StepBwd,test="LRT")
anova(m_sature,m_StepFwd,test="LRT")Pour les 2 modèles m_StepBwdFwd, m_StepBwd: on rejette de manière très nette l’hypothèse H0. Et on accepte la nullité des paramètres du modèle : m_sature qui ne sont pas dans les modèles step. On préfèrera donc m_StepBwdFwd, m_StepBwd. Par contre pour le modèle m_StepFwd le test n’est pas significatif. En effet au seuil de 6% on préfèrera le modèle saturé.
4.5.2 Les modèles obtenus à partir de l’analyse PCA
Le modèle issue de l’analyse PCA
# Modèle obtenue avec l'analyse PCA - Avec l'option quali.sup=[1:4]
m_PCA<-glm(pluie.demain ~ 1 +
WGust_mean + WGust_max #1er axe
+ WS.10m_mean + WS.80m_mean + WS.900mb_mean #1er axe
+ WS.10m_max + WS.80m_max + WS.900mb_max #1er axe
+ Tot.CC_mean + Humidity_mean + Humidity_min + #2eme axe
+ WS.10m_min + WS.80m_min + #2eme axe
+ Temp_mean + Temp_max + Temp_min + #3eme axe
+ WDir.80m_mean + WDir.10m_mean + WDir.900mb_mean + #4eme axe
+ Tot.CC_min + CC.low_min + CC.mid_min #5eme axe
+ CC.mid_max + MSL_max + MSL_mean + MSL_min #divers
+ Sunshine_sum + Tot.Prcp_sum + Snow_sum + Radiation_sum, #divers
family = binomial, data = meteoTrain )
summary(m_PCA)##
## Call:
## glm(formula = pluie.demain ~ 1 + WGust_mean + WGust_max + WS.10m_mean +
## WS.80m_mean + WS.900mb_mean + WS.10m_max + WS.80m_max + WS.900mb_max +
## Tot.CC_mean + Humidity_mean + Humidity_min + +WS.10m_min +
## WS.80m_min + +Temp_mean + Temp_max + Temp_min + +WDir.80m_mean +
## WDir.10m_mean + WDir.900mb_mean + +Tot.CC_min + CC.low_min +
## CC.mid_min + CC.mid_max + MSL_max + MSL_mean + MSL_min +
## Sunshine_sum + Tot.Prcp_sum + Snow_sum + Radiation_sum, family = binomial,
## data = meteoTrain)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.4848 -0.8890 -0.2178 0.9211 2.3324
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 8.483e+01 1.174e+01 7.227 4.94e-13 ***
## WGust_mean 3.122e-02 2.951e-02 1.058 0.29011
## WGust_max 1.224e-02 1.485e-02 0.825 0.40961
## WS.10m_mean -5.978e-02 9.352e-02 -0.639 0.52267
## WS.80m_mean -9.866e-02 6.442e-02 -1.532 0.12560
## WS.900mb_mean -3.314e-03 1.519e-02 -0.218 0.82735
## WS.10m_max 2.713e-02 3.238e-02 0.838 0.40217
## WS.80m_max 1.185e-03 2.674e-02 0.044 0.96464
## WS.900mb_max 7.873e-03 9.487e-03 0.830 0.40663
## Tot.CC_mean 2.749e-03 7.318e-03 0.376 0.70719
## Humidity_mean -1.029e-02 1.792e-02 -0.574 0.56564
## Humidity_min 1.473e-02 1.411e-02 1.044 0.29642
## WS.10m_min 1.074e-01 5.448e-02 1.971 0.04872 *
## WS.80m_min -6.750e-04 3.874e-02 -0.017 0.98610
## Temp_mean 1.883e-01 1.546e-01 1.218 0.22311
## Temp_max 2.224e-02 9.098e-02 0.244 0.80686
## Temp_min -2.232e-01 8.101e-02 -2.755 0.00587 **
## WDir.80m_mean -8.179e-03 5.200e-03 -1.573 0.11570
## WDir.10m_mean 6.514e-03 5.016e-03 1.299 0.19402
## WDir.900mb_mean 6.206e-03 1.350e-03 4.596 4.30e-06 ***
## Tot.CC_min 2.009e-02 8.295e-03 2.422 0.01543 *
## CC.low_min -2.200e-02 8.402e-03 -2.619 0.00882 **
## CC.mid_min 6.779e-03 1.059e-02 0.640 0.52218
## CC.mid_max 9.541e-03 2.158e-03 4.421 9.81e-06 ***
## MSL_max -1.647e-01 7.249e-02 -2.272 0.02307 *
## MSL_mean 2.033e-01 1.280e-01 1.588 0.11227
## MSL_min -1.246e-01 6.569e-02 -1.898 0.05776 .
## Sunshine_sum -4.711e-04 8.269e-04 -0.570 0.56884
## Tot.Prcp_sum 2.971e-02 2.246e-02 1.323 0.18580
## Snow_sum 8.477e-01 4.795e-01 1.768 0.07710 .
## Radiation_sum 2.922e-05 8.693e-05 0.336 0.73680
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1719.0 on 1239 degrees of freedom
## Residual deviance: 1354.7 on 1209 degrees of freedom
## AIC: 1416.7
##
## Number of Fisher Scoring iterations: 6Les 2 modèles obtenus par minimsation du critère AIC de la méthode step à partir de m_PCA
(m_StepBwdFwd_PCA et m_StepBwd_PCA sont identiques)
#Modèle step backward/forward: m_StepBwdFwd_PCA à partir de m_PCA
m_StepBwdFwd_PCA<-glm(formula = pluie.demain ~ 1 + WGust_max + WS.80m_mean + WS.10m_min + Temp_mean +
Temp_min + WDir.900mb_mean + Tot.CC_min + CC.low_min + CC.mid_max +
MSL_max + MSL_mean + MSL_min + Sunshine_sum + Tot.Prcp_sum + Snow_sum,
family = binomial, data = meteoTrain)
summary(m_StepBwdFwd_PCA)##
## Call:
## glm(formula = pluie.demain ~ 1 + WGust_max + WS.80m_mean + WS.10m_min +
## Temp_mean + Temp_min + WDir.900mb_mean + Tot.CC_min + CC.low_min +
## CC.mid_max + MSL_max + MSL_mean + MSL_min + Sunshine_sum +
## Tot.Prcp_sum + Snow_sum, family = binomial, data = meteoTrain)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.4704 -0.8986 -0.2133 0.9206 2.2978
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 84.3069781 11.2135349 7.518 5.55e-14 ***
## WGust_max 0.0303372 0.0086109 3.523 0.000427 ***
## WS.80m_mean -0.0881162 0.0221422 -3.980 6.90e-05 ***
## WS.10m_min 0.0953245 0.0295656 3.224 0.001263 **
## Temp_mean 0.2220405 0.0510736 4.347 1.38e-05 ***
## Temp_min -0.2290999 0.0538483 -4.255 2.09e-05 ***
## WDir.900mb_mean 0.0049273 0.0010664 4.621 3.83e-06 ***
## Tot.CC_min 0.0225040 0.0075161 2.994 0.002753 **
## CC.low_min -0.0236141 0.0081172 -2.909 0.003624 **
## CC.mid_max 0.0096377 0.0020870 4.618 3.87e-06 ***
## MSL_max -0.1684431 0.0711393 -2.368 0.017895 *
## MSL_mean 0.2206233 0.1249573 1.766 0.077465 .
## MSL_min -0.1374704 0.0634936 -2.165 0.030380 *
## Sunshine_sum -0.0007125 0.0003833 -1.859 0.063076 .
## Tot.Prcp_sum 0.0358339 0.0215180 1.665 0.095853 .
## Snow_sum 0.7675947 0.4462083 1.720 0.085385 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1719.0 on 1239 degrees of freedom
## Residual deviance: 1362.2 on 1224 degrees of freedom
## AIC: 1394.2
##
## Number of Fisher Scoring iterations: 6#Modèle step forward: m_StepFwd_PCA à partir de m_PCA
m_StepFwd_PCA<-glm(formula = pluie.demain ~ MSL_min + CC.mid_max + WDir.900mb_mean +
Snow_sum + CC.mid_min + MSL_max + MSL_mean + WS.900mb_max +
WS.80m_mean + WS.10m_min + WGust_max + WDir.80m_mean,
family = binomial, data = meteoTrain)
summary(m_StepFwd_PCA)##
## Call:
## glm(formula = pluie.demain ~ MSL_min + CC.mid_max + WDir.900mb_mean +
## Snow_sum + CC.mid_min + MSL_max + MSL_mean + WS.900mb_max +
## WS.80m_mean + WS.10m_min + WGust_max + WDir.80m_mean, family = binomial,
## data = meteoTrain)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.2865 -0.9121 -0.2727 0.9528 2.3796
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 91.542919 10.746683 8.518 < 2e-16 ***
## MSL_min -0.109850 0.061527 -1.785 0.0742 .
## CC.mid_max 0.011204 0.001851 6.052 1.43e-09 ***
## WDir.900mb_mean 0.004855 0.001145 4.241 2.23e-05 ***
## Snow_sum 0.830187 0.444741 1.867 0.0619 .
## CC.mid_min 0.016493 0.009101 1.812 0.0700 .
## MSL_max -0.139374 0.068238 -2.042 0.0411 *
## MSL_mean 0.157736 0.120555 1.308 0.1907
## WS.900mb_max 0.007511 0.005862 1.281 0.2000
## WS.80m_mean -0.089385 0.021896 -4.082 4.46e-05 ***
## WS.10m_min 0.079620 0.028617 2.782 0.0054 **
## WGust_max 0.025399 0.010154 2.501 0.0124 *
## WDir.80m_mean -0.001961 0.001338 -1.466 0.1426
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1719.0 on 1239 degrees of freedom
## Residual deviance: 1390.8 on 1227 degrees of freedom
## AIC: 1416.8
##
## Number of Fisher Scoring iterations: 6Comparaison avec le modèle “saturé”:
anova(m_sature,m_StepBwdFwd_PCA,test="LRT")
anova(m_sature,m_StepFwd_PCA,test="LRT")On conservera m_StepBwdFwd_PCA. Mais ce n’est pas le cas pour m_StepFwd_PCA.
Comparaison avec le modèle “PCA”:
anova(m_PCA,m_StepBwdFwd_PCA,test="LRT")
anova(m_PCA,m_StepFwd_PCA,test="LRT")Là encore m_StepBwdFwd_PCA est meilleur que m_PCA mais ce n’est pas le cas de m_StepFwd_PCA. On va donc éliminer m_StepFwd_PCA
Les 2 modèles obtenus à partir de la méthode leaps par minimsation du critère R2 ajusté puis Cp de Mallows
(obtenue à partir de la préselection par PCA) m_StepBwdFwd_PCA et m_Cp_PCA sont identiques on a donc éliminé m_Cp_PCA: Cp de Mallows
#Modèle leaps critère R2 ajusté: m_R2aj_PCA à partir de m_PCA
m_R2aj_PCA<-glm(pluie.demain ~
Temp_mean + MSL_mean + Tot.Prcp_sum + Snow_sum + Sunshine_sum +
WDir.10m_mean + WS.80m_mean + WDir.80m_mean + WDir.900mb_mean +
Temp_min + MSL_max + MSL_min + Tot.CC_min + CC.mid_max + CC.low_min +
WS.10m_min + WS.900mb_max + WGust_max, family = binomial, data = meteoTrain )
summary(m_R2aj_PCA)##
## Call:
## glm(formula = pluie.demain ~ Temp_mean + MSL_mean + Tot.Prcp_sum +
## Snow_sum + Sunshine_sum + WDir.10m_mean + WS.80m_mean + WDir.80m_mean +
## WDir.900mb_mean + Temp_min + MSL_max + MSL_min + Tot.CC_min +
## CC.mid_max + CC.low_min + WS.10m_min + WS.900mb_max + WGust_max,
## family = binomial, data = meteoTrain)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.4856 -0.8923 -0.2156 0.9126 2.3530
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 85.0898498 11.2695300 7.550 4.34e-14 ***
## Temp_mean 0.2152700 0.0527172 4.083 4.44e-05 ***
## MSL_mean 0.2048491 0.1257596 1.629 0.10334
## Tot.Prcp_sum 0.0332619 0.0216302 1.538 0.12411
## Snow_sum 0.8148584 0.4534238 1.797 0.07232 .
## Sunshine_sum -0.0006931 0.0003855 -1.798 0.07220 .
## WDir.10m_mean 0.0067393 0.0049196 1.370 0.17072
## WS.80m_mean -0.0969541 0.0227306 -4.265 2.00e-05 ***
## WDir.80m_mean -0.0078760 0.0051575 -1.527 0.12674
## WDir.900mb_mean 0.0058789 0.0012877 4.565 4.99e-06 ***
## Temp_min -0.2209433 0.0557242 -3.965 7.34e-05 ***
## MSL_max -0.1653974 0.0713801 -2.317 0.02050 *
## MSL_min -0.1254476 0.0640655 -1.958 0.05022 .
## Tot.CC_min 0.0230780 0.0075770 3.046 0.00232 **
## CC.mid_max 0.0096849 0.0020906 4.633 3.61e-06 ***
## CC.low_min -0.0236804 0.0081811 -2.895 0.00380 **
## WS.10m_min 0.0971550 0.0296524 3.276 0.00105 **
## WS.900mb_max 0.0071674 0.0060624 1.182 0.23710
## WGust_max 0.0254466 0.0102911 2.473 0.01341 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1719.0 on 1239 degrees of freedom
## Residual deviance: 1358.3 on 1221 degrees of freedom
## AIC: 1396.3
##
## Number of Fisher Scoring iterations: 6#Modèle leaps critère Cp Mallows: m_Cp_PCA à partir de m_PCA identique à m_StepBwdFwd_PCA
m_Cp_PCA<-glm(pluie.demain ~
Temp_mean + MSL_mean + Tot.Prcp_sum + Sunshine_sum + WS.80m_mean + WDir.900mb_mean
+ Temp_min + MSL_max + MSL_min + Tot.CC_min + CC.mid_max + CC.low_min
+ WGust_max + WS.10m_min ,family = binomial, data = meteoTrain )
summary(m_Cp_PCA)##
## Call:
## glm(formula = pluie.demain ~ Temp_mean + MSL_mean + Tot.Prcp_sum +
## Sunshine_sum + WS.80m_mean + WDir.900mb_mean + Temp_min +
## MSL_max + MSL_min + Tot.CC_min + CC.mid_max + CC.low_min +
## WGust_max + WS.10m_min, family = binomial, data = meteoTrain)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.4877 -0.9006 -0.2210 0.9196 2.2715
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 84.8014601 11.1555733 7.602 2.92e-14 ***
## Temp_mean 0.2197111 0.0510527 4.304 1.68e-05 ***
## MSL_mean 0.2350596 0.1234275 1.904 0.056854 .
## Tot.Prcp_sum 0.0374880 0.0215099 1.743 0.081364 .
## Sunshine_sum -0.0006897 0.0003827 -1.802 0.071508 .
## WS.80m_mean -0.0844408 0.0218956 -3.857 0.000115 ***
## WDir.900mb_mean 0.0050376 0.0010629 4.739 2.14e-06 ***
## Temp_min -0.2318763 0.0538680 -4.305 1.67e-05 ***
## MSL_max -0.1768833 0.0701809 -2.520 0.011722 *
## MSL_min -0.1439403 0.0629903 -2.285 0.022306 *
## Tot.CC_min 0.0230426 0.0075220 3.063 0.002189 **
## CC.mid_max 0.0097456 0.0020838 4.677 2.91e-06 ***
## CC.low_min -0.0234188 0.0080864 -2.896 0.003779 **
## WGust_max 0.0297735 0.0085777 3.471 0.000518 ***
## WS.10m_min 0.0903657 0.0292203 3.093 0.001984 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1719.0 on 1239 degrees of freedom
## Residual deviance: 1366.7 on 1225 degrees of freedom
## AIC: 1396.7
##
## Number of Fisher Scoring iterations: 4Comparaison avec le modèle “saturé” d’origine:
anova(m_sature,m_R2aj_PCA,test="LRT")
anova(m_PCA,m_R2aj_PCA,test="LRT")
anova(m_sature,m_Cp_PCA,test="LRT")
anova(m_PCA,m_Cp_PCA,test="LRT")On préfèrera la aussi m_R2aj_PCA et m_Cp_PCA au modèle saturé: m_sature On conserve le modèle m_R2aj_PCA et m_Cp_PCA.
4.5.3 Comparaison des différents modèles avec les critères BIC et AIC
Tous ces modèles ont un meilleur BIC et AIC que le modèle saturé. A la vue des critères BIC et AIC les modèles (hors modèles initiaux: saturé et PCA) semblent très proches. Excepté pour le modèle : m_StepFwd qui est en retrait par rapport aux autres.
| Modèle | AIC | BIC |
|---|---|---|
| m_sature | 1431.849 | 1662.378 |
| m_StepBwdFwd | 1392.372 | 1474.338 |
| m_StepBwd | 1393.164 | 1470.007 |
| m_StepFwd | 1413.588 | 1480.185 |
| m_PCA | 1416.738 | 1575.547 |
| m_StepBwdFwd_PCA | 1394.18 | 1476.146 |
| m_R2aj_PCA | 1396.347 | 1493.681 |
| m_Cp_PCA | 1396.679 | 1473.522 |
On va essayer de les départager en utilsant la validation croisée.
5. Validation croisé - Performance des modèles
Nos données de test sont incomplètes, il manque le résultat pluie.demain. On ne peut donc pas effectuer une validation croisée àpartir de ces données. Pour mener à bien cette comparaison il faut séparer (à nouveaux) les données d’entrainement (meteoTrain) en 2 jeux de données : dataTest et dataTrain. Puis recalibrer les modèles obtenus à partir du nouveau jeux de données: dataTrain. Et effectuer les tests sur le jeux dataTest. Ceci n’est qu’une approximation et donc à considérer avec retenue. En effet les nouvelles données d’entrainement seraient succeptibles de produire de nouveaux modèles.
5.0 Initialisation
On doit recharger tous les modèles pour les calibrer avec les nouvelles données. On effectue la prédiction pour chacun des modèles, à partir de la fonction : predict appliquée aux données de tests (dataTest) Avec l’option : type=“reponse” pour bien avoir la probabilité de prédiction.
5.1 Probabilité estimé d’avoir pluie.demain=TRUE
## Sature StepBwdFwd StepBwd StepFwd PCA StepBwdFwd_PCA R2aj_PCA Cp_PCA
## 2 0.516 0.551 0.502 0.492 0.417 0.436 0.415 0.430
## 3 0.267 0.337 0.295 0.279 0.253 0.308 0.277 0.298
## 6 0.819 0.743 0.765 0.785 0.854 0.808 0.800 0.803
## 9 0.520 0.732 0.731 0.722 0.510 0.690 0.666 0.699
## 14 0.253 0.314 0.318 0.289 0.280 0.256 0.263 0.244
## 18 0.238 0.264 0.298 0.362 0.258 0.311 0.306 0.298
## 19 0.218 0.221 0.195 0.202 0.194 0.184 0.184 0.1775.2 Estimation au seuil de 0.5
On confronte les probabilités obtenues en 5.1 aux seuil de 0.5 Dés que la prévision dépasse 50% on prédit qu’il pleuvra demain.
pred_0.5 <- apply(pred_proba >=0.5, 2, factor,labels=c("FALSE","TRUE"))
head(pred_0.5,n=7)## Sature StepBwdFwd StepBwd StepFwd PCA StepBwdFwd_PCA R2aj_PCA
## 2 "TRUE" "TRUE" "TRUE" "FALSE" "FALSE" "FALSE" "FALSE"
## 3 "FALSE" "FALSE" "FALSE" "FALSE" "FALSE" "FALSE" "FALSE"
## 6 "TRUE" "TRUE" "TRUE" "TRUE" "TRUE" "TRUE" "TRUE"
## 9 "TRUE" "TRUE" "TRUE" "TRUE" "TRUE" "TRUE" "TRUE"
## 14 "FALSE" "FALSE" "FALSE" "FALSE" "FALSE" "FALSE" "FALSE"
## 18 "FALSE" "FALSE" "FALSE" "FALSE" "FALSE" "FALSE" "FALSE"
## 19 "FALSE" "FALSE" "FALSE" "FALSE" "FALSE" "FALSE" "FALSE"
## Cp_PCA
## 2 "FALSE"
## 3 "FALSE"
## 6 "TRUE"
## 9 "TRUE"
## 14 "FALSE"
## 18 "FALSE"
## 19 "FALSE"5.3 Pourcentage de réussite des modèles par rapport à l’observation
On compare maintenant les valeurs prédites aux valeurs observées. On fait mieux que la moyenne, par contre on a très peu gagné par rapport au modèle saturé, entre 1% et 2%. La méthode employée s’est révélée assez peu performante.
## Sature StepBwdFwd StepBwd StepFwd PCA StepBwdFwd_PCA R2aj_PCA
## 1 70.2479 73.1405 73.1405 73.5537 71.9008 72.7273 73.5537
## Cp_PCA Observe
## 1 73.5537 1005.4 Courbes ROC
Le meilleur estiamteur aura une aire sous la courbe le plus proche possible de 1. La courbe idéale serait perpendiculaire aux abscisse du point origine jusqu’au point (0.1) puis prallèle jusqu’au point (1,1) Si bien que l’aire sous cette courbe serait égale à 1.
6. Choix et validation du modèle
Le modèle retenu est le modèle: m_StepBwdFwd. Obtenu à partir du modèle complet (m_sature) par méthode pas à pas progessive en utilsant la fonction step du package MASS qui minimise le critère AIC. De tous nos modèles c’est celui qui a le plus peit AIC. Et l’un des critère BIC les plus faible. La validation croisé est assez satisfaisante, avec un taux de succès de l’ordre de 70%. Mais la grande majorité des modèles construits sont très proches et il est difficile de les départager. Il faudrait pouvoir réaliser une véritable validation croisée et utiliser différentes estimation de l’écart aux données observées (Normes L1, L2..).
6.1 Corrélations entre variables explicatives
On remarque qu’un nombre important de variable est peu corrélée. Ce qui est une très bonne chose. L’information apportées par les vqriables du modèle est bien complémentaire.
Par contre nos graphiques ne permettent pas de détecter la multicollinéarité. Pour celà on peut utiliser la fonction VIF (variance inflation factors) du package Faraway: VIF = 1/(1−R2) Plus le VIF est grand, plus la variable est corrélée avec les autres. VIF ∈ [1, +∞[
require(faraway)
vif(m_StepBwdFwd)## Temp_mean Snow_sum CC.low_mean WS.80m_mean
## 175.622139 88.271852 14.440561 38.737045
## WDir.900mb_mean Temp_min MSL_max Tot.CC_max
## 7.740385 163.168211 8.910557 11.766211
## Tot.CC_min CC.high_max CC.mid_max CC.low_min
## 37.159657 9.872790 14.713964 27.890780
## WS.10m_min WGust_max WS.900mb_max
## 19.077444 29.358666 20.873122Problème de colinéarité entre les variables Temp_mean et Temp_min. Le résultat générale semble assez moyen. L’utilisation du R2, et donc de la fonction VIF elle même est possiblement mal adapté au modèles glm.
6.2 Résidus et résidus partiels
Tous les résultats sont basés sur des hypothèses fondamentales liés au terme d’erreur qui résume les informations absentes du modèle. Hypothèses liées au terme d’erreur :
- Les erreurs sont centrées E(erreur) = 0.
- La distribution est gaussienne
- La variance est constante (σ2)
- Les erreurs (Erreur_i)1≤i≤n sont indépendantes.
Pour vérifier ces hypothèses on va regarder des graphiques des résidus qui sont des approximations de ces erreurs.
- L’hypothèse d’indépendance est considérée comme vérifiée lorsque chaque donnée correspond à un échantillonnage indépendant ou à une expérience physique menée dans des conditions indépendantes. Ce qui est le cas ici.
- On peut souvent se passer de l’hypothèse de normalité général quand le nombre le données est important. Ce qui est bien le cas.
- La plupart des graphiques montrent bien une répartition linéaire centrèe en 0.
6.3 Points influents / abérrants
Il conviendrait d’analyser les points influents/levier et déterminer pourquoi ils le sont. Ici on se contente juste de les visualiser.
which(influence(m_StepBwdFwd)$hat>0.15)## 717 1033
## 713 1029
which(cooks.distance(m_StepBwdFwd)>0.03)## named integer(0)7. Prévision - Modèle glm - Tests
7.1 Chargement des données de test
meteoTest= read.csv("meteo.test.csv")
meteoTest<-renameCol(meteoTest)
summary(meteoTest)## X Year Month Day Hour
## Min. : 221 Min. :2011 Min. : 1.000 Min. : 1.0 Min. :0
## 1st Qu.: 888 1st Qu.:2012 1st Qu.: 4.000 1st Qu.: 8.0 1st Qu.:0
## Median :1619 Median :2014 Median : 7.000 Median :15.0 Median :0
## Mean :1603 Mean :2014 Mean : 6.665 Mean :15.5 Mean :0
## 3rd Qu.:2331 3rd Qu.:2016 3rd Qu.:10.000 3rd Qu.:23.0 3rd Qu.:0
## Max. :2939 Max. :2018 Max. :12.000 Max. :31.0 Max. :0
## Minute Temp_mean Humidity_mean MSL_mean
## Min. :0 Min. :-1.750 Min. :44.00 Min. : 994
## 1st Qu.:0 1st Qu.: 6.045 1st Qu.:63.88 1st Qu.:1012
## Median :0 Median :11.880 Median :71.58 Median :1017
## Mean :0 Mean :11.797 Mean :70.96 Mean :1017
## 3rd Qu.:0 3rd Qu.:16.730 3rd Qu.:77.56 3rd Qu.:1023
## Max. :0 Max. :29.960 Max. :93.50 Max. :1047
## Tot.Prcp_sum Snow_sum Tot.CC_mean CC.high_mean
## Min. : 0.000 Min. :0.00000 Min. : 0.00 Min. : 0.00
## 1st Qu.: 0.000 1st Qu.:0.00000 1st Qu.: 22.85 1st Qu.: 4.02
## Median : 0.000 Median :0.00000 Median : 57.83 Median :15.50
## Mean : 1.803 Mean :0.03022 Mean : 52.49 Mean :22.08
## 3rd Qu.: 2.150 3rd Qu.:0.00000 3rd Qu.: 78.47 3rd Qu.:35.56
## Max. :21.100 Max. :3.50000 Max. :100.00 Max. :92.54
## CC.mid_mean CC.low_mean Sunshine_sum Radiation_sum
## Min. : 0.000 Min. : 0.000 Min. : 0.0 Min. : 284.8
## 1st Qu.: 3.855 1st Qu.: 9.685 1st Qu.:115.2 1st Qu.:1966.0
## Median : 28.500 Median : 40.120 Median :352.1 Median :3259.2
## Mean : 33.284 Mean : 40.974 Mean :365.3 Mean :3792.6
## 3rd Qu.: 55.315 3rd Qu.: 66.020 3rd Qu.:584.8 3rd Qu.:5368.9
## Max. :100.000 Max. :100.000 Max. :947.8 Max. :8370.5
## WS.10m_mean WDir.10m_mean WS.80m_mean WDir.80m_mean
## Min. : 3.130 Min. : 60.88 Min. : 3.540 Min. : 59.99
## 1st Qu.: 6.345 1st Qu.:160.94 1st Qu.: 8.705 1st Qu.:171.19
## Median : 8.990 Median :211.71 Median :12.110 Median :219.48
## Mean :10.742 Mean :204.34 Mean :14.306 Mean :210.76
## 3rd Qu.:13.255 3rd Qu.:248.13 3rd Qu.:18.450 3rd Qu.:255.60
## Max. :37.480 Max. :328.74 Max. :48.220 Max. :323.77
## WS.900mb_mean WDir.900mb_mean WGust_mean Temp_max
## Min. : 4.00 Min. : 58.55 Min. : 3.110 Min. : 1.080
## 1st Qu.:13.03 1st Qu.:144.74 1st Qu.: 9.655 1st Qu.: 9.615
## Median :19.68 Median :235.44 Median :13.340 Median :16.220
## Mean :24.44 Mean :206.73 Mean :16.685 Mean :15.974
## 3rd Qu.:30.48 3rd Qu.:265.50 3rd Qu.:20.730 3rd Qu.:21.120
## Max. :85.57 Max. :315.95 Max. :55.900 Max. :36.700
## Temp_min Humidity_max Humidity_min MSL_max
## Min. :-5.460 Min. : 65.00 Min. :20.00 Min. : 996.5
## 1st Qu.: 2.420 1st Qu.: 82.00 1st Qu.:44.00 1st Qu.:1014.9
## Median : 7.810 Median : 88.00 Median :54.00 Median :1019.1
## Mean : 7.784 Mean : 87.11 Mean :53.85 Mean :1020.0
## 3rd Qu.:12.455 3rd Qu.: 94.00 3rd Qu.:62.00 3rd Qu.:1025.2
## Max. :24.290 Max. :100.00 Max. :89.00 Max. :1049.9
## MSL_min Tot.CC_max Tot.CC_min CC.high_max
## Min. : 991.7 Min. : 0.00 Min. : 0.000 Min. : 0.00
## 1st Qu.:1008.7 1st Qu.:100.00 1st Qu.: 0.000 1st Qu.: 24.00
## Median :1014.6 Median :100.00 Median : 0.000 Median :100.00
## Mean :1014.4 Mean : 89.87 Mean : 9.618 Mean : 66.11
## 3rd Qu.:1019.6 3rd Qu.:100.00 3rd Qu.: 3.000 3rd Qu.:100.00
## Max. :1043.6 Max. :100.00 Max. :100.000 Max. :100.00
## CC.high_min CC.mid_max CC.mid_min CC.low_max
## Min. : 0.0000 Min. : 0.0 Min. : 0.000 Min. : 0.00
## 1st Qu.: 0.0000 1st Qu.: 55.5 1st Qu.: 0.000 1st Qu.:100.00
## Median : 0.0000 Median :100.0 Median : 0.000 Median :100.00
## Mean : 0.7445 Mean : 76.1 Mean : 2.357 Mean : 80.95
## 3rd Qu.: 0.0000 3rd Qu.:100.0 3rd Qu.: 0.000 3rd Qu.:100.00
## Max. :27.0000 Max. :100.0 Max. :100.000 Max. :100.00
## CC.low_min WS.10m_max WS.10m_min WS.80m_max
## Min. : 0.00 Min. : 5.24 Min. : 0.000 Min. : 5.35
## 1st Qu.: 0.00 1st Qu.:11.61 1st Qu.: 0.800 1st Qu.:17.43
## Median : 0.00 Median :17.21 Median : 1.940 Median :23.21
## Mean : 4.96 Mean :18.82 Mean : 3.549 Mean :24.81
## 3rd Qu.: 0.00 3rd Qu.:22.32 3rd Qu.: 4.800 3rd Qu.:29.44
## Max. :100.00 Max. :52.34 Max. :21.420 Max. :66.49
## WS.80m_min WS.900mb_max WS.900mb_min WGust_max
## Min. : 0.000 Min. : 7.59 Min. : 0.360 Min. : 6.48
## 1st Qu.: 1.140 1st Qu.: 24.36 1st Qu.: 2.585 1st Qu.:18.36
## Median : 2.310 Median : 34.52 Median : 7.240 Median :25.20
## Mean : 5.059 Mean : 40.35 Mean :11.122 Mean :28.84
## 3rd Qu.: 6.840 3rd Qu.: 50.92 3rd Qu.:14.050 3rd Qu.:34.02
## Max. :30.300 Max. :109.60 Max. :61.350 Max. :83.16
## WGust_min
## Min. : 0.000
## 1st Qu.: 2.160
## Median : 4.320
## Mean : 6.888
## 3rd Qu.: 8.640
## Max. :39.2407.2 prévision
on a la proba p(x) = P(Y=Oui/X=x)
prev_m_StepBwdFwd <- predict(m_StepBwdFwd,newdata = meteoTest, type = "response" )
summary(prev_m_StepBwdFwd)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.01102 0.29638 0.51445 0.49633 0.70214 0.98007prev_m_StepBwdFwd[1:20]## 1 2 3 4 5 6
## 0.90957013 0.33140101 0.04412599 0.41587243 0.38560748 0.28939204
## 7 8 9 10 11 12
## 0.13397579 0.55712220 0.66902023 0.39742119 0.70504950 0.30677207
## 13 14 15 16 17 18
## 0.40668011 0.60971473 0.76466549 0.75684454 0.88144080 0.66608398
## 19 20
## 0.60523983 0.57159316Si on décide de prédire qu’il pleuvra demain avec un seuil de 0.70. Alors les jours: 1, 11, 15, 16, 17, seront prédit avec de la pluie. Si on abaisse le seuil à 0.6 on aura en plus les jours: 9, 14, 18, 19 qui seront prédit avec de la pluie.
Prévision pour les autres modèles:
prev_sature <- predict(m_sature,newdata = meteoTest, type = "response" )
summary(prev_sature)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.01227 0.28745 0.48520 0.49218 0.69081 0.97538prev_m_StepBwd <- predict(m_StepBwd,newdata = meteoTest, type = "response" )
summary(prev_m_StepBwd)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.009114 0.299541 0.497602 0.497221 0.703840 0.981449prev_m_StepFwd <- predict(m_StepFwd,newdata = meteoTest, type = "response" )
summary(prev_m_StepFwd)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.02167 0.30381 0.50007 0.49691 0.69598 0.98094prev_PCA <- predict(m_PCA,newdata = meteoTest, type = "response" )
summary(prev_PCA)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.01334 0.28486 0.48751 0.48769 0.68076 0.97590prev_m_StepBwdFwd_PCA <- predict(m_StepBwdFwd_PCA,newdata = meteoTest, type = "response" )
summary(prev_m_StepBwdFwd_PCA)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.01178 0.29224 0.48526 0.49103 0.70718 0.96437prev_m_R2aj_PCA <- predict(m_R2aj_PCA,newdata = meteoTest, type = "response" )
summary(prev_m_R2aj_PCA)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.01187 0.29852 0.48335 0.48765 0.68104 0.96765prev_m_Cp_PCA <- predict(m_Cp_PCA,newdata = meteoTest, type = "response" )
summary(prev_m_Cp_PCA)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.01218 0.30019 0.48726 0.49129 0.69167 0.963667.3 Probabilité estimé d’avoir pluie.demain=TRUE
## prev_sature prev_StepBwdFwd prev_StepBwd prev_StepFwd prev_PCA
## 1 0.952 0.910 0.893 0.910 0.916
## 2 0.365 0.331 0.359 0.547 0.386
## 3 0.063 0.044 0.050 0.052 0.074
## 4 0.401 0.416 0.375 0.302 0.324
## 5 0.274 0.386 0.382 0.336 0.319
## 6 0.249 0.289 0.269 0.257 0.266
## 7 0.146 0.134 0.094 0.183 0.123
## prev_StepBwdFwd.PCA prev_R2aj.PCA prev_CpMallow.PCA
## 1 0.883 0.893 0.894
## 2 0.432 0.397 0.445
## 3 0.048 0.074 0.050
## 4 0.333 0.309 0.347
## 5 0.380 0.388 0.377
## 6 0.291 0.280 0.299
## 7 0.096 0.099 0.1017.4 Estimation au seuil de 0.5
On confronte les probabilités obtenues en 6.1 aux seuil de 0.5 Dés que la prévision dépasse 50% on prédit qu’il pleuvra demain.
## prev_sature prev_StepBwdFwd prev_StepBwd prev_StepFwd prev_PCA
## 1 "TRUE" "TRUE" "TRUE" "TRUE" "TRUE"
## 2 "FALSE" "FALSE" "FALSE" "TRUE" "FALSE"
## 3 "FALSE" "FALSE" "FALSE" "FALSE" "FALSE"
## 4 "FALSE" "FALSE" "FALSE" "FALSE" "FALSE"
## 5 "FALSE" "FALSE" "FALSE" "FALSE" "FALSE"
## 6 "FALSE" "FALSE" "FALSE" "FALSE" "FALSE"
## 7 "FALSE" "FALSE" "FALSE" "FALSE" "FALSE"
## prev_StepBwdFwd.PCA prev_R2aj.PCA prev_CpMallow.PCA
## 1 "TRUE" "TRUE" "TRUE"
## 2 "FALSE" "FALSE" "FALSE"
## 3 "FALSE" "FALSE" "FALSE"
## 4 "FALSE" "FALSE" "FALSE"
## 5 "FALSE" "FALSE" "FALSE"
## 6 "FALSE" "FALSE" "FALSE"
## 7 "FALSE" "FALSE" "FALSE"7.5 Concaténation aux données météoTest et constitution du fichier de sortie
meteoTestOut<-meteoTest
meteoTestOut<-renameCol_inv(meteoTestOut)
res<-cbind(meteoTestOut,m_pred_0.5)
res_All<-cbind(meteoTestOut,m_pred_0.5_All)
head(res_All,n=7)
write.csv(res,"C:\\MeteoPrevision.csv", quote=FALSE)
write.csv(res_All,"C:\\MeteoPrevision_All.csv", quote=FALSE)Conclusion
On retient le modèle: m_StepBwdFwd. C’est le modèle qui a le meilleur AIC et l’un des critère BIC les plus faible. Sinon la pluspart des modèles sont très proches en terme de critère AIC et BIC. Seuls les modèles de type step forward sont un peu en retrait. Les modèles m_StepBwdFwd, m_StepBwd et m_StepBwdFwd_PCA ont des critères BIC et AIC très voisins.
Il semble normal, d’obtenir un moins bon critère AIC pour notre modèle: m_StepBwdFwd_PCA que pour le modèle : m_StepBwdFwd. En effet ces deux modèles sont obtenus par minimisation du critère AIC à partir de la fonction step du package MASS. Ils utilisent donc la même méthode. La différence vient de la pré-selection PCA, par rapport au modèle saturé. Ainsi notre pré-selection “manuelle” s’est avérée moins performante. Cependant l’écart est très faible.
Tous ces modèles construit à partir du modèle saturé et du modèle issue de l’analyse PCA: m_PCA ont de meilleurs critères que les modèles d’origine. Notamment de meilleurs test anova. Pour les comparaisons anova et BIc, AIC on pourra revenir au §4.5.
Les librairies leaps (avec le critère R2 ajusté et Cp de Mallow - cf 4.3 ) et glmulti (cf. 4.4) n’ont pas apporté l’amélioration attendue. Vraissemblablement la pré-sélection effectuée à partir de l’analyse PCA était mal adaptée. La libraire glmulti et sa méthode associée est difficilement utilsable, de part le coût de calcul. Le maximun de 30 variables est difficilement utilisable, il engendre plusieurs heures de calcul avec 2^29 test de modèles.
On note la compléxité de l’analyse météorologique, avec un très grand nombre de vraiables et la difficulté à réduire leur nombre.
La mise en place d’une (certaine) validation croisée nous a permis d’avoir une certaine idée du pourcentage de réussite. Qui se situe autour de 68%-70%. Nous n’avons pas pu utiliser nos données de test pour une validation croisé et une évaluation de performance. Nos données de test étant incomplètes: il manque le résultat pluie.demain. Pour mener à bien cette comparaison (§5) il a fallu séparer (à nouveaux) les données d’entrainement en 2 jeux de données : dataTest et dataTrain. Recalibrer les modèles obtenus à partir du nouveau jeux de données: dataTrain. Puis effectuer les tests sur le jeux dataTest. Ceci n’est qu’une approximation et donc à considérer avec retenue. En effet les nouvelles données d’entrainement sont succeptible de produire de nouveaux modèles.
Le fichier résultat : MeteoPrevision_All.csv contient les données de test avec des colonnes supplémentaires qui sont les prédictions de tous ces modèles. Le ficher résultat ; MeteoPrevision.csv contient les données de test avec la colonne pluie.demain qui est la prédiction à partir du modèle retenu: m_StepBkwdFwd.